Calibration

The calibration package provides tools for fitting financial models to market data, estimating parameters from historical data, and calculating implied volatility surfaces.

BSMIVSolver

BSMIVSolver(max_iter: int = 20, tolerance: float = 1e-06, *, bisection_iter: int = 60, iv_bounds: tuple[float, float] = (0.0001, 5.0))

High-performance, vectorized implied-volatility solver.

Safeguarded Newton-Raphson (vega floor, step clipping, bracket tracking) on every option at once, followed by bisection on the ones Newton did not converge. Prices outside the no-arbitrage band (below the discounted intrinsic value or above the forward/strike bound) get NaN.

Initializes the BSM implied volatility solver.

Parameters:

Name	Type	Description	Default
max_iter	int	Newton iterations, by default 20.	20
tolerance	float	Absolute price tolerance for convergence, by default 1e-6.	1e-06
bisection_iter	int	Bisection iterations for options Newton could not solve, default 60.	60
iv_bounds	tuple[float, float]	Search interval for the volatility, by default (1e-4, 5.0).	(0.0001, 5.0)

Source code in src/optpricing/calibration/vectorized_bsm_iv.py

def __init__(
    self,
    max_iter: int = 20,
    tolerance: float = 1e-6,
    *,
    bisection_iter: int = 60,
    iv_bounds: tuple[float, float] = (1e-4, 5.0),
):
    """
    Initializes the BSM implied volatility solver.

    Parameters
    ----------
    max_iter : int, optional
        Newton iterations, by default 20.
    tolerance : float, optional
        Absolute price tolerance for convergence, by default 1e-6.
    bisection_iter : int, optional
        Bisection iterations for options Newton could not solve, default 60.
    iv_bounds : tuple[float, float], optional
        Search interval for the volatility, by default (1e-4, 5.0).
    """
    self.max_iter = max_iter
    self.tolerance = tolerance
    self.bisection_iter = bisection_iter
    self.iv_bounds = iv_bounds

solve

solve(target_prices: ndarray, options: DataFrame, stock: Stock, rate: Rate, *, initial_guess: ndarray | None = None) -> np.ndarray

Calculates implied volatility for an array of options.

Parameters:

Name	Type	Description	Default
target_prices	ndarray	An array of market prices for which to find the implied volatility.	required
options	DataFrame	A DataFrame of option contracts, must include 'strike', 'maturity', and 'optionType' columns.	required
stock	Stock	The underlying asset's properties.	required
rate	Rate	The risk-free rate structure.	required
initial_guess	ndarray | None	Starting volatilities (e.g. yesterday's or the market's implied vols); defaults to 0.2 everywhere.	None

Returns:

Type	Description
ndarray	Implied volatilities aligned with target_prices; NaN where the price admits no Black-Scholes volatility.

Source code in src/optpricing/calibration/vectorized_bsm_iv.py

def solve(
    self,
    target_prices: np.ndarray,
    options: pd.DataFrame,
    stock: Stock,
    rate: Rate,
    *,
    initial_guess: np.ndarray | None = None,
) -> np.ndarray:
    """
    Calculates implied volatility for an array of options.

    Parameters
    ----------
    target_prices : np.ndarray
        An array of market prices for which to find the implied volatility.
    options : pd.DataFrame
        A DataFrame of option contracts, must include 'strike', 'maturity',
        and 'optionType' columns.
    stock : Stock
        The underlying asset's properties.
    rate : Rate
        The risk-free rate structure.
    initial_guess : np.ndarray | None, optional
        Starting volatilities (e.g. yesterday's or the market's implied
        vols); defaults to 0.2 everywhere.

    Returns
    -------
    np.ndarray
        Implied volatilities aligned with ``target_prices``; ``NaN`` where
        the price admits no Black-Scholes volatility.
    """
    target = np.asarray(target_prices, dtype=float)
    S, q = stock.spot, stock.dividend
    K = options["strike"].to_numpy(dtype=float)
    T = options["maturity"].to_numpy(dtype=float)
    if callable(rate.rate):
        r = np.asarray([rate.get_rate(t) for t in T], dtype=float)
    else:
        r = np.full(T.shape, float(rate.get_rate()))
    is_call = options["optionType"].to_numpy() == "call"
    lo, hi = self.iv_bounds

    # No-arbitrage band: intrinsic (discounted) <= price <= S e^{-qT} or K e^{-rT}
    fwd_leg, strike_leg = S * np.exp(-q * T), K * np.exp(-r * T)
    intrinsic = np.where(is_call, fwd_leg - strike_leg, strike_leg - fwd_leg)
    upper = np.where(is_call, fwd_leg, strike_leg)
    solvable = (target > np.maximum(intrinsic, 0.0)) & (target < upper)
    solvable &= np.isfinite(target) & (T > 0) & (K > 0)

    iv = (
        np.full(target.shape, 0.2)
        if initial_guess is None
        else np.array(initial_guess, dtype=float)
    )
    iv = np.where(np.isfinite(iv), np.clip(iv, lo, hi), 0.2)
    active = solvable.copy()
    low_b, high_b = np.full(target.shape, lo), np.full(target.shape, hi)

    for _ in range(self.max_iter):
        if not active.any():
            break
        prices, vega = bsm_prices_and_vegas(iv, S, K, T, r, q, is_call)
        error = prices - target
        converged = np.abs(error) < self.tolerance
        active &= ~converged
        # Price is increasing in vol: tighten the bracket on the way.
        low_b = np.where(active & (error < 0), np.maximum(low_b, iv), low_b)
        high_b = np.where(active & (error > 0), np.minimum(high_b, iv), high_b)
        step = error / np.maximum(vega, 1e-12)
        new_iv = iv - np.clip(step, -1.0, 1.0)
        inside = (new_iv > low_b) & (new_iv < high_b)
        new_iv = np.where(inside, new_iv, 0.5 * (low_b + high_b))
        iv = np.where(active, new_iv, iv)

    # Bisection for whatever Newton left over (deep ITM/OTM, tiny vegas).
    if active.any():
        idx = np.flatnonzero(active)
        a, b = low_b[idx], high_b[idx]
        for _ in range(self.bisection_iter):
            mid = 0.5 * (a + b)
            prices, _ = bsm_prices_and_vegas(
                mid, S, K[idx], T[idx], r[idx], q, is_call[idx]
            )
            too_high = prices > target[idx]
            a, b = np.where(too_high, a, mid), np.where(too_high, mid, b)
        iv[idx] = 0.5 * (a + b)
        prices, _ = bsm_prices_and_vegas(
            iv[idx], S, K[idx], T[idx], r[idx], q, is_call[idx]
        )
        solvable[idx] &= np.abs(prices - target[idx]) < 1e-3 * np.maximum(
            target[idx], 1.0
        )

    return np.where(solvable, iv, np.nan)

CalibrationResult dataclass

CalibrationResult(params: dict[str, float], loss: float, rmse_iv: float, rmse_price: float, n_options: int, n_evaluations: int, n_starts: int, runtime: float, success: bool, message: str = '', starts: list[dict[str, Any]] = list())

Outcome of :meth:Calibrator.fit.

Attributes:

Name	Type	Description
params	dict[str, float]	Calibrated (and frozen) parameters.
loss	float	Final weighted mean squared residual (IV² or price² units).
rmse_iv	float	Unweighted RMSE of the implied-volatility residuals at the optimum (NaN when no option admits an implied volatility).
rmse_price	float	Unweighted RMSE of the price residuals at the optimum.
n_options	int	Size of the calibration set.
n_evaluations	int	Objective evaluations (all starts).
n_starts	int	Starting points tried.
runtime	float	Wall-clock seconds.
success	bool	Whether at least one local optimisation converged.
message	str	Optimiser message of the best run.
starts	list[dict[str, Any]]	Per-start records (params, loss, kind).

CalibrationSettings dataclass

CalibrationSettings(liquidity: LiquidityFilter = LiquidityFilter(), objective: str = 'iv', weighting: str = 'spread', n_starts: int = 4, warm_start: bool = True, historical_prior: bool = True, feller_penalty: float = 0.0, time_budget: float | None = None, name: str = 'default')

Everything a calibration run needs besides the data and the model.

Parameters:

Name	Type	Description	Default
liquidity	LiquidityFilter	Selection of the calibration set.	LiquidityFilter()
objective	('iv', 'price')	Residuals in Black-Scholes implied-volatility or in price space.	"iv"
weighting	('spread', 'vega', 'none')	Residual weights: inverse bid-ask spread, inverse vega (price space only; approximates the IV objective) or equal.	"spread"
n_starts	int	Number of starting points (Latin hypercube in the bounds) besides the recipe's initial guess and the warm start; 1 = local optimisation only.	4
warm_start	bool	Seed the optimiser with the previous calibration (when one is given).	True
historical_prior	bool	Use historical-return estimates for jump parameters as a start point.	True
feller_penalty	float	Soft penalty weight on vol_of_vol^2 - 2 kappa theta for Heston-type models (0 disables it).	0.0
time_budget	float | None	Wall-clock seconds per model (None = unlimited).	None
name	str	Label used in tuning reports.	'default'

replace

replace(**changes: Any) -> CalibrationSettings

A copy with some fields changed (liquidity may be a dict of filter fields to update).

Source code in src/optpricing/calibration/market_filters.py

def replace(self, **changes: Any) -> CalibrationSettings:
    """A copy with some fields changed (``liquidity`` may be a dict of
    filter fields to update)."""
    liq = changes.pop("liquidity", None)
    data = self.to_dict()
    data.update(changes)
    if isinstance(liq, LiquidityFilter):
        data["liquidity"] = liq.to_dict()
    elif isinstance(liq, dict):
        data["liquidity"] = {**data["liquidity"], **liq}
    return CalibrationSettings.from_dict(data)

Calibrator

Calibrator(model: BaseModel, market_data: DataFrame, stock: Stock, rate: Rate, *, objective: str = 'iv', weighting: str = 'spread', allow_mc: bool = False, feller_penalty: float = 0.0, mc_paths: int = 20000, mc_steps: int = 50, random_state: int = 0)

Fits a model's parameters to an option chain.

The residuals are taken in Black-Scholes implied-volatility space by default (objective="iv"), which puts every strike and maturity on the same footing; objective="price" uses price differences. Residuals can be weighted by the inverse bid-ask spread (default), the inverse vega (price space) or equally. Parameters are optimised in a unit cube mapped onto the bounds with L-BFGS-B, from several starting points (initial guess, optional warm start and Latin-hypercube samples); invalid parameter sets (model ValueError) and optional Feller violations are penalised.

Initializes the Calibrator.

Parameters:

Name	Type	Description	Default
model	BaseModel	The financial model to be calibrated (e.g., HestonModel).	required
market_data	DataFrame	Option quotes: 'strike', 'maturity', 'optionType', 'marketPrice' and optionally 'bid'/'ask' (for spread weights) and 'impliedVolatility' (warm start for the IV solver).	required
stock	Stock	The underlying asset's properties.	required
rate	Rate	The risk-free rate structure.	required
objective	('iv', 'price')	Residual space, by default "iv".	"iv"
weighting	('spread', 'vega', 'none')	Residual weights, by default "spread".	"spread"
allow_mc	bool	Permit Monte Carlo pricing for simulation-only models.	False
feller_penalty	float	Weight of the soft Feller penalty (Heston/Bates), by default 0.	0.0
mc_paths	int	Monte Carlo settings when allow_mc is used.	20000
mc_steps	int	Monte Carlo settings when allow_mc is used.	20000
random_state	int	Seed for the Latin-hypercube starting points, by default 0.	0

Raises:

Type	Description
NotImplementedError	If the model cannot be priced on a chain (no closed form, no characteristic function and allow_mc is False).

Source code in src/optpricing/calibration/calibrator.py

def __init__(
    self,
    model: BaseModel,
    market_data: pd.DataFrame,
    stock: Stock,
    rate: Rate,
    *,
    objective: str = "iv",
    weighting: str = "spread",
    allow_mc: bool = False,
    feller_penalty: float = 0.0,
    mc_paths: int = 20_000,
    mc_steps: int = 50,
    random_state: int = 0,
):
    """
    Initializes the Calibrator.

    Parameters
    ----------
    model : BaseModel
        The financial model to be calibrated (e.g., HestonModel).
    market_data : pd.DataFrame
        Option quotes: 'strike', 'maturity', 'optionType', 'marketPrice'
        and optionally 'bid'/'ask' (for spread weights) and
        'impliedVolatility' (warm start for the IV solver).
    stock : Stock
        The underlying asset's properties.
    rate : Rate
        The risk-free rate structure.
    objective : {"iv", "price"}, optional
        Residual space, by default "iv".
    weighting : {"spread", "vega", "none"}, optional
        Residual weights, by default "spread".
    allow_mc : bool, optional
        Permit Monte Carlo pricing for simulation-only models.
    feller_penalty : float, optional
        Weight of the soft Feller penalty (Heston/Bates), by default 0.
    mc_paths, mc_steps : int, optional
        Monte Carlo settings when ``allow_mc`` is used.
    random_state : int, optional
        Seed for the Latin-hypercube starting points, by default 0.

    Raises
    ------
    NotImplementedError
        If the model cannot be priced on a chain (no closed form, no
        characteristic function and ``allow_mc`` is False).
    """
    if objective not in ("iv", "price"):
        raise ValueError("objective must be 'iv' or 'price'.")
    if weighting not in ("spread", "vega", "none"):
        raise ValueError("weighting must be 'spread', 'vega' or 'none'.")
    if not supports_vectorized_pricing(model, allow_mc=allow_mc):
        raise NotImplementedError(
            f"{model.name} cannot be calibrated: no closed form or "
            "characteristic function (pass allow_mc=True to use Monte Carlo)."
        )
    self.model = model
    self.stock = stock
    self.rate = rate
    self.objective = objective
    self.weighting = weighting
    self.allow_mc = allow_mc
    self.feller_penalty = float(feller_penalty)
    self.mc_kwargs = {"mc_paths": mc_paths, "mc_steps": mc_steps}
    self.random_state = random_state
    self.technique = select_fastest_technique(model)
    self.iv_solver = BSMIVSolver()
    self.market_data = self._prepare(market_data)
    self.result: CalibrationResult | None = None
    self.n_evaluations = 0
    logger.info(
        "Calibrator: %s on %d options, objective=%s, weighting=%s",
        model.name,
        len(self.market_data),
        objective,
        weighting,
    )

fit

fit(initial_guess: dict[str, float], bounds: dict[str, tuple], frozen_params: dict[str, float] | None = None, *, n_starts: int = 1, warm_start: dict[str, float] | None = None, n_local: int = 3, time_budget: float | None = None, maxiter: int = 200, log_scale: Iterable[str] = ()) -> dict[str, float]

Performs the calibration.

Parameters:

Name	Type	Description	Default
initial_guess	dict[str, float]	Initial guesses for the parameters to be fitted.	required
bounds	dict[str, tuple]	(min, max) bounds per fitted parameter; missing bounds default to ±50 % around the guess.	required
frozen_params	dict[str, float] | None	Parameters held constant. Defaults to None.	None
n_starts	int	Total starting points: the guess, the warm start (if given) and Latin-hypercube samples filling up to n_starts. Default 1.	1
warm_start	dict[str, float] | None	Previous calibration used as an extra start point.	None
n_local	int	How many of the best-scoring starts get a full local optimisation.	3
time_budget	float | None	Wall-clock seconds; further starts are skipped once exceeded.	None
maxiter	int	L-BFGS-B iteration cap per start, by default 200.	200
log_scale	Iterable[str]	Parameters mapped logarithmically onto their (positive) bounds, for scales spanning orders of magnitude (e.g. the CEV sigma, whose natural size is vol * S^(1-gamma)).	()

Returns:

Type	Description
dict[str, float]	The full set of calibrated and frozen parameters. Details are in self.result (:class:CalibrationResult).

Source code in src/optpricing/calibration/calibrator.py

def fit(
    self,
    initial_guess: dict[str, float],
    bounds: dict[str, tuple],
    frozen_params: dict[str, float] | None = None,
    *,
    n_starts: int = 1,
    warm_start: dict[str, float] | None = None,
    n_local: int = 3,
    time_budget: float | None = None,
    maxiter: int = 200,
    log_scale: Iterable[str] = (),
) -> dict[str, float]:
    """
    Performs the calibration.

    Parameters
    ----------
    initial_guess : dict[str, float]
        Initial guesses for the parameters to be fitted.
    bounds : dict[str, tuple]
        (min, max) bounds per fitted parameter; missing bounds default to
        ±50 % around the guess.
    frozen_params : dict[str, float] | None, optional
        Parameters held constant. Defaults to None.
    n_starts : int, optional
        Total starting points: the guess, the warm start (if given) and
        Latin-hypercube samples filling up to ``n_starts``. Default 1.
    warm_start : dict[str, float] | None, optional
        Previous calibration used as an extra start point.
    n_local : int, optional
        How many of the best-scoring starts get a full local optimisation.
    time_budget : float | None, optional
        Wall-clock seconds; further starts are skipped once exceeded.
    maxiter : int, optional
        L-BFGS-B iteration cap per start, by default 200.
    log_scale : Iterable[str], optional
        Parameters mapped logarithmically onto their (positive) bounds, for
        scales spanning orders of magnitude (e.g. the CEV ``sigma``, whose
        natural size is ``vol * S^(1-gamma)``).

    Returns
    -------
    dict[str, float]
        The full set of calibrated and frozen parameters. Details are in
        ``self.result`` (:class:`CalibrationResult`).
    """
    t0 = time.perf_counter()
    self.n_evaluations = 0
    frozen = dict(frozen_params or {})
    names = [p for p in initial_guess if p not in frozen]
    if not names:
        self.result = CalibrationResult(
            params=frozen,
            loss=self.loss(frozen) if len(self.market_data) else float("nan"),
            rmse_iv=float("nan"),
            rmse_price=float("nan"),
            n_options=len(self.market_data),
            n_evaluations=self.n_evaluations,
            n_starts=0,
            runtime=time.perf_counter() - t0,
            success=True,
            message="nothing to fit",
        )
        return frozen
    if len(self.market_data) == 0:
        raise ValueError("No options to calibrate to.")

    lo, hi = self._bounds_arrays(names, initial_guess, bounds)
    is_log = np.array([p in set(log_scale) for p in names])
    if np.any(is_log & (lo <= 0)):
        raise ValueError("log_scale parameters need strictly positive bounds.")
    log_lo = np.log(np.where(is_log, lo, 1.0))
    log_hi = np.log(np.where(is_log, hi, np.e))

    def to_unit(x: Any) -> np.ndarray:
        x = np.asarray(x, dtype=float)
        linear = (x - lo) / (hi - lo)
        with np.errstate(divide="ignore", invalid="ignore"):
            logarithmic = (np.log(np.where(is_log, x, 1.0)) - log_lo) / (
                log_hi - log_lo
            )
        return np.where(is_log, logarithmic, linear)

    def from_unit(u: np.ndarray) -> np.ndarray:
        u = np.clip(u, 0.0, 1.0)
        linear = lo + u * (hi - lo)
        logarithmic = np.exp(log_lo + u * (log_hi - log_lo))
        return np.where(is_log, logarithmic, linear)

    def objective(u: np.ndarray) -> float:
        return self.loss({**frozen, **dict(zip(names, from_unit(u)))})

    # -- starting points
    starts: list[tuple[np.ndarray, str]] = [
        (np.clip(to_unit([initial_guess[p] for p in names]), 0, 1), "guess")
    ]
    if warm_start and all(p in warm_start for p in names):
        starts.append(
            (np.clip(to_unit([warm_start[p] for p in names]), 0, 1), "warm")
        )
    n_random = max(n_starts - len(starts), 0)
    if n_random:
        sampler = qmc.LatinHypercube(d=len(names), seed=self.random_state)
        for u in sampler.random(n_random):
            starts.append((u, "lhs"))

    scored = []
    for u, kind in starts:
        scored.append((objective(u), u, kind))
    scored.sort(key=lambda s: s[0])

    # -- local optimisation from the best starts
    records: list[dict[str, Any]] = []
    best_u, best_loss, best_msg, success = scored[0][1], scored[0][0], "", False
    for i, (loss0, u0, kind) in enumerate(scored[: max(n_local, 1)]):
        if time_budget is not None and i > 0:
            if time.perf_counter() - t0 > time_budget:
                logger.info("Time budget exhausted after %d start(s).", i)
                break
        res = minimize(
            objective,
            x0=u0,
            method="L-BFGS-B",
            bounds=[(0.0, 1.0)] * len(names),
            options={"maxiter": maxiter},
        )
        records.append(
            {
                "kind": kind,
                "start_loss": float(loss0),
                "loss": float(res.fun),
                "params": dict(zip(names, from_unit(res.x))),
                "success": bool(res.success),
            }
        )
        if res.fun < best_loss:
            best_u, best_loss, best_msg = res.x, float(res.fun), str(res.message)
        success = success or bool(res.success)

    final = {**frozen, **dict(zip(names, from_unit(best_u)))}
    rmse_iv, rmse_price = self._final_rmse(final)
    self.result = CalibrationResult(
        params=final,
        loss=best_loss,
        rmse_iv=rmse_iv,
        rmse_price=rmse_price,
        n_options=len(self.market_data),
        n_evaluations=self.n_evaluations,
        n_starts=len(starts),
        runtime=time.perf_counter() - t0,
        success=success,
        message=best_msg,
        starts=records,
    )
    logger.info(
        "%s calibrated: loss=%.3e rmse_iv=%.4f rmse_price=%.4f "
        "(%d evals, %d starts, %.1fs)",
        self.model.name,
        best_loss,
        rmse_iv,
        rmse_price,
        self.n_evaluations,
        len(starts),
        self.result.runtime,
    )
    return final

loss

loss(params: dict[str, float]) -> float

Weighted mean squared residual (plus penalties) for a parameter set.

Source code in src/optpricing/calibration/calibrator.py

def loss(self, params: dict[str, float]) -> float:
    """Weighted mean squared residual (plus penalties) for a parameter set."""
    self.n_evaluations += 1
    try:
        model = self.model.with_params(**params)
    except ValueError as e:
        logger.debug("Invalid params %s (%s)", params, e)
        return self._invalid_loss
    try:
        res = self._residuals(model)
    except (ValueError, FloatingPointError, ZeroDivisionError) as e:
        logger.debug("Pricing failed for %s (%s)", params, e)
        return self._invalid_loss
    if not np.all(np.isfinite(res)):
        return self._invalid_loss
    value = float(np.sum(self._weights * res * res) / np.sum(self._weights))
    return value + self._penalty(params)

LiquidityFilter dataclass

LiquidityFilter(min_maturity_days: float = 7.0, max_maturity_years: float = 2.0, moneyness: tuple[float, float] = (0.7, 1.3), max_relative_spread: float = 0.25, min_open_interest: int = 1, require_bid: bool = True, otm_only: bool = True, iv_bounds: tuple[float, float] = (0.02, 3.0))

Rules selecting the options a model is fitted to (or evaluated on).

Parameters:

Name	Type	Description	Default
min_maturity_days	float	Maturity window; default 7 days to 2 years.	7.0
max_maturity_years	float	Maturity window; default 7 days to 2 years.	7.0
moneyness	tuple[float, float]	Allowed forward moneyness K / F with F = S e^{(r-q)T}; default (0.7, 1.3).	(0.7, 1.3)
max_relative_spread	float	Maximum (ask - bid) / mid; default 0.25. Ignored without quotes.	0.25
min_open_interest	int	Minimum open interest (or volume when open interest is missing); 1.	1
require_bid	bool	Drop quotes with a zero bid; default True.	True
otm_only	bool	Keep only out-of-the-money options (puts below the forward, calls above), the most liquid side of the chain; default True.	True
iv_bounds	tuple[float, float]	Market implied volatilities outside this range are treated as data errors; default (0.02, 3.0).	(0.02, 3.0)

evaluation

evaluation() -> LiquidityFilter

The wider companion filter used to evaluate a fit: moneyness widened by 0.1 on each side and twice the spread tolerance, same maturity window, liquidity floor and OTM rule (in-the-money American puts carry an early-exercise premium no European model should be judged on).

Source code in src/optpricing/calibration/market_filters.py

def evaluation(self) -> LiquidityFilter:
    """
    The wider companion filter used to *evaluate* a fit: moneyness widened
    by 0.1 on each side and twice the spread tolerance, same maturity
    window, liquidity floor and OTM rule (in-the-money American puts carry
    an early-exercise premium no European model should be judged on).
    """
    lo, hi = self.moneyness
    return LiquidityFilter(
        min_maturity_days=self.min_maturity_days,
        max_maturity_years=self.max_maturity_years,
        moneyness=(max(lo - 0.1, 0.05), hi + 0.1),
        max_relative_spread=min(2.0 * self.max_relative_spread, 1.0),
        min_open_interest=self.min_open_interest,
        require_bid=self.require_bid,
        otm_only=self.otm_only,
        iv_bounds=self.iv_bounds,
    )

to_dict

to_dict() -> dict[str, Any]

JSON-ready dict (tuples become lists, as JSON would store them).

Source code in src/optpricing/calibration/market_filters.py

def to_dict(self) -> dict[str, Any]:
    """JSON-ready dict (tuples become lists, as JSON would store them)."""
    data = asdict(self)
    for key in ("moneyness", "iv_bounds"):
        data[key] = list(data[key])
    return data

VectorizedIntegrationIVSolver

VectorizedIntegrationIVSolver(max_iter: int = 20, tolerance: float = 1e-07, upper_bound: float = 1000.0)

A high-performance, vectorized Secant method solver for implied volatility for any model that supports a characteristic function.

Initializes the vectorized IV solver.

Parameters:

Name	Type	Description	Default
max_iter	int	Maximum number of iterations for the Secant method, by default 20.	20
tolerance	float	Error tolerance for convergence, by default 1e-7.	1e-07
upper_bound	float	The upper limit for the numerical integration, by default 1000.0.	1000.0

Source code in src/optpricing/calibration/vectorized_integration_iv.py

def __init__(
    self,
    max_iter: int = 20,
    tolerance: float = 1e-7,
    upper_bound: float = 1000.0,
):
    """
    Initializes the vectorized IV solver.

    Parameters
    ----------
    max_iter : int, optional
        Maximum number of iterations for the Secant method, by default 20.
    tolerance : float, optional
        Error tolerance for convergence, by default 1e-7.
    upper_bound : float, optional
        The upper limit for the numerical integration, by default 1000.0.
    """
    self.max_iter = max_iter
    self.tolerance = tolerance
    self.upper_bound = upper_bound

solve

solve(target_prices: ndarray, options: DataFrame, model: BaseModel, rate: Rate) -> np.ndarray

Calculates implied volatility for an array of options and prices.

This method uses a vectorized Secant root-finding algorithm. The pricing at each step is performed using a vectorized version of the Gil-Pelaez inversion formula.

Parameters:

Name	Type	Description	Default
target_prices	ndarray	An array of market prices for which to find the implied volatility.	required
options	DataFrame	A DataFrame of option contracts.	required
model	BaseModel	The financial model to use for pricing.	required
rate	Rate	The risk-free rate structure.	required

Returns:

Type	Description
ndarray	An array of calculated implied volatilities.

Source code in src/optpricing/calibration/vectorized_integration_iv.py

def solve(
    self,
    target_prices: np.ndarray,
    options: pd.DataFrame,
    model: BaseModel,
    rate: Rate,
) -> np.ndarray:
    """
    Calculates implied volatility for an array of options and prices.

    This method uses a vectorized Secant root-finding algorithm. The pricing
    at each step is performed using a vectorized version of the Gil-Pelaez
    inversion formula.

    Parameters
    ----------
    target_prices : np.ndarray
        An array of market prices for which to find the implied volatility.
    options : pd.DataFrame
        A DataFrame of option contracts.
    model : BaseModel
        The financial model to use for pricing.
    rate : Rate
        The risk-free rate structure.

    Returns
    -------
    np.ndarray
        An array of calculated implied volatilities.
    """
    iv0 = np.full_like(target_prices, 0.20)
    iv1 = np.full_like(target_prices, 0.25)

    p0 = self._price_vectorized(iv0, options, model, rate)
    p1 = self._price_vectorized(iv1, options, model, rate)

    f0 = p0 - target_prices
    f1 = p1 - target_prices

    for _ in range(self.max_iter):
        if np.all(np.abs(f1) < self.tolerance):
            break
        denom = f1 - f0
        denom[np.abs(denom) < 1e-12] = 1e-12
        iv_next = iv1 - f1 * (iv1 - iv0) / denom
        iv_next = np.clip(iv_next, 1e-4, 5.0)
        iv0, iv1 = iv1, iv_next
        f0 = f1
        p1 = self._price_vectorized(iv1, options, model, rate)
        f1 = p1 - target_prices

    return iv1

VolatilitySurface

VolatilitySurface(option_data: DataFrame)

A class to compute and hold market and model-implied volatility surfaces.

This class takes a DataFrame of option market data and provides methods to calculate the Black-Scholes implied volatility (IV) for each option, either from its market price or from a price generated by a financial model.

Initializes the VolatilitySurface.

Parameters:

Name	Type	Description	Default
option_data	DataFrame	A DataFrame containing market prices of options. Must include 'strike', 'maturity', 'marketPrice', 'optionType', and 'expiry' columns.	required

Raises:

Type	Description
ValueError	If any of the required columns are missing from option_data.

Source code in src/optpricing/calibration/iv_surface.py

def __init__(self, option_data: pd.DataFrame):
    """
    Initializes the VolatilitySurface.

    Parameters
    ----------
    option_data : pd.DataFrame
        A DataFrame containing market prices of options. Must include
        'strike', 'maturity', 'marketPrice', 'optionType', and 'expiry' columns.

    Raises
    ------
    ValueError
        If any of the required columns are missing from `option_data`.
    """
    required_cols = ["strike", "maturity", "marketPrice", "optionType", "expiry"]
    if not all(col in option_data.columns for col in required_cols):
        msg = (
            "Input option_data is missing one "
            f"of the required columns: {required_cols}"
        )
        raise ValueError(msg)

    self.data = option_data[required_cols].copy()
    self.surface: pd.DataFrame | None = None
    self.iv_solver = BSMIVSolver()

calculate_market_iv

calculate_market_iv(stock: Stock, rate: Rate) -> VolatilitySurface

Calculates the market implied volatility surface from market prices.

This method inverts the Black-Scholes formula for each option's market price to find the corresponding implied volatility. The results are stored in the self.surface DataFrame.

Parameters:

Name	Type	Description	Default
stock	Stock	The underlying asset's properties.	required
rate	Rate	The risk-free rate structure.	required

Returns:

Type	Description
VolatilitySurface	The same instance of the class, allowing for method chaining.

Source code in src/optpricing/calibration/iv_surface.py

def calculate_market_iv(
    self,
    stock: Stock,
    rate: Rate,
) -> VolatilitySurface:
    """
    Calculates the market implied volatility surface from market prices.

    This method inverts the Black-Scholes formula for each option's market
    price to find the corresponding implied volatility. The results are
    stored in the `self.surface` DataFrame.

    Parameters
    ----------
    stock : Stock
        The underlying asset's properties.
    rate : Rate
        The risk-free rate structure.

    Returns
    -------
    VolatilitySurface
        The same instance of the class, allowing for method chaining.
    """
    logger.info("Calculating market implied volatility surface...")
    market_ivs = self._calculate_ivs(stock, rate, self.data["marketPrice"])
    self.surface = self.data.copy()
    self.surface["iv"] = market_ivs
    self.surface = self.surface.dropna()
    self.surface = self.surface[
        (self.surface["iv"] > 1e-4) & (self.surface["iv"] < 2.0)
    ]
    return self

calculate_model_iv

calculate_model_iv(stock: Stock, rate: Rate, model: BaseModel, technique: BaseTechnique = None, *, allow_mc: bool = False) -> VolatilitySurface

Calculates a model's implied volatility surface.

This method first prices every option in the dataset using the provided model (vectorised closed form, characteristic function or, with allow_mc, Monte Carlo). It then inverts the Black-Scholes formula for each of these model prices to generate the model-implied volatility surface.

Parameters:

Name	Type	Description	Default
stock	Stock	The underlying asset's properties.	required
rate	Rate	The risk-free rate structure.	required
model	BaseModel	The financial model to generate prices from.	required
technique	BaseTechnique	Unused (kept for backwards compatibility); the chain pricer picks the fastest exact method itself.	None
allow_mc	bool	Permit Monte Carlo pricing for simulation-only models.	False

Returns:

Type	Description
VolatilitySurface	The same instance of the class, allowing for method chaining.

Source code in src/optpricing/calibration/iv_surface.py

def calculate_model_iv(
    self,
    stock: Stock,
    rate: Rate,
    model: BaseModel,
    technique: BaseTechnique = None,
    *,
    allow_mc: bool = False,
) -> VolatilitySurface:
    """
    Calculates a model's implied volatility surface.

    This method first prices every option in the dataset using the provided
    model (vectorised closed form, characteristic function or, with
    ``allow_mc``, Monte Carlo). It then inverts the Black-Scholes formula
    for each of these model prices to generate the model-implied
    volatility surface.

    Parameters
    ----------
    stock : Stock
        The underlying asset's properties.
    rate : Rate
        The risk-free rate structure.
    model : BaseModel
        The financial model to generate prices from.
    technique : BaseTechnique
        Unused (kept for backwards compatibility); the chain pricer picks
        the fastest exact method itself.
    allow_mc : bool, optional
        Permit Monte Carlo pricing for simulation-only models.

    Returns
    -------
    VolatilitySurface
        The same instance of the class, allowing for method chaining.
    """
    logger.info("Calculating %s implied volatility surface...", model.name)

    model_prices = price_options_vectorized(
        options_df=self.data,
        stock=stock,
        model=model,
        rate=rate,
        allow_mc=allow_mc,
        **model.params,
    )

    model_ivs = self._calculate_ivs(
        stock, rate, pd.Series(model_prices, index=self.data.index)
    )
    self.surface = self.data.copy()
    self.surface["iv"] = model_ivs
    self.surface = self.surface.dropna()
    return self

annotate_chain

annotate_chain(df: DataFrame, stock: Stock, rate: Rate, *, iv_solver: BSMIVSolver | None = None) -> pd.DataFrame

Adds pricing helpers to a chain: forward, moneyness (K/F), log_moneyness, relative_spread and iv_market (Black-Scholes implied volatility of the mid, NaN where none exists).

Parameters:

Name	Type	Description	Default
df	DataFrame	Options with strike, maturity, optionType, marketPrice and optionally bid/ask.	required
stock	(Stock, Rate)	Market parameters used for the forward and the implied volatilities.	required
rate	(Stock, Rate)	Market parameters used for the forward and the implied volatilities.	required
iv_solver	BSMIVSolver | None	Solver instance (default settings if omitted).	None

Returns:

Type	Description
DataFrame	A copy of df with the extra columns and a fresh integer index.

Source code in src/optpricing/calibration/market_filters.py

def annotate_chain(
    df: pd.DataFrame,
    stock: Stock,
    rate: Rate,
    *,
    iv_solver: BSMIVSolver | None = None,
) -> pd.DataFrame:
    """
    Adds pricing helpers to a chain: ``forward``, ``moneyness`` (``K/F``),
    ``log_moneyness``, ``relative_spread`` and ``iv_market`` (Black-Scholes
    implied volatility of the mid, NaN where none exists).

    Parameters
    ----------
    df : pd.DataFrame
        Options with ``strike``, ``maturity``, ``optionType``, ``marketPrice``
        and optionally ``bid``/``ask``.
    stock, rate : Stock, Rate
        Market parameters used for the forward and the implied volatilities.
    iv_solver : BSMIVSolver | None, optional
        Solver instance (default settings if omitted).

    Returns
    -------
    pd.DataFrame
        A copy of ``df`` with the extra columns and a fresh integer index.
    """
    out = df.reset_index(drop=True).copy()
    T = out["maturity"].to_numpy(dtype=float)
    if callable(rate.rate):
        r = np.asarray([rate.get_rate(t) for t in T], dtype=float)
    else:
        r = np.full(T.shape, float(rate.get_rate()))
    forward = stock.spot * np.exp((r - stock.dividend) * T)
    out["forward"] = forward
    out["moneyness"] = out["strike"].to_numpy(dtype=float) / forward
    out["log_moneyness"] = np.log(out["moneyness"].to_numpy())
    if {"bid", "ask"} <= set(out.columns):
        mid = out["marketPrice"].to_numpy(dtype=float)
        with np.errstate(divide="ignore", invalid="ignore"):
            out["relative_spread"] = (out["ask"] - out["bid"]).to_numpy() / mid
    else:
        out["relative_spread"] = np.nan
    solver = iv_solver or BSMIVSolver()
    guess = (
        out["impliedVolatility"].to_numpy(dtype=float)
        if "impliedVolatility" in out.columns
        else None
    )
    out["iv_market"] = solver.solve(
        out["marketPrice"].to_numpy(dtype=float), out, stock, rate, initial_guess=guess
    )
    return out

bucket_metrics

bucket_metrics(table: DataFrame) -> pd.DataFrame

RMSE of the IV and price errors by maturity and by moneyness bucket.

Returns:

Type	Description
DataFrame	Columns bucket_type ("maturity" / "moneyness"), bucket, n_options, rmse_iv, rmse_price.

Source code in src/optpricing/calibration/metrics.py

def bucket_metrics(table: pd.DataFrame) -> pd.DataFrame:
    """
    RMSE of the IV and price errors by maturity and by moneyness bucket.

    Returns
    -------
    pd.DataFrame
        Columns ``bucket_type`` ("maturity" / "moneyness"), ``bucket``,
        ``n_options``, ``rmse_iv``, ``rmse_price``.
    """
    if table.empty:
        return pd.DataFrame(
            columns=["bucket_type", "bucket", "n_options", "rmse_iv", "rmse_price"]
        )
    rows = []
    groups = {
        "maturity": maturity_bucket(table["maturity"]),
        "moneyness": moneyness_bucket(table["moneyness"]),
    }
    for kind, labels in groups.items():
        for label, grp in table.groupby(labels, observed=True):
            rows.append(
                {
                    "bucket_type": kind,
                    "bucket": str(label),
                    "n_options": int(len(grp)),
                    "rmse_iv": _rmse(grp["iv_error"]),
                    "rmse_price": _rmse(grp["price_error"]),
                }
            )
    return pd.DataFrame(rows)

evaluate_model

evaluate_model(model: BaseModel, data: DataFrame, stock: Stock, rate: Rate, *, allow_mc: bool = False, **pricing_kwargs: Any) -> tuple[dict[str, float], pd.DataFrame, pd.DataFrame]

One-call evaluation: scalar metrics, bucket table and residual table.

Parameters:

Name	Type	Description	Default
model	BaseModel	See :func:residual_table.	required
data	BaseModel	See :func:residual_table.	required
stock	BaseModel	See :func:residual_table.	required
rate	BaseModel	See :func:residual_table.	required
allow_mc	bool	Permit the Monte Carlo pricing route.	False

Returns:

Type	Description
tuple[dict[str, float], DataFrame, DataFrame]	(summary, buckets, residuals).

Source code in src/optpricing/calibration/metrics.py

def evaluate_model(
    model: BaseModel,
    data: pd.DataFrame,
    stock: Stock,
    rate: Rate,
    *,
    allow_mc: bool = False,
    **pricing_kwargs: Any,
) -> tuple[dict[str, float], pd.DataFrame, pd.DataFrame]:
    """
    One-call evaluation: scalar metrics, bucket table and residual table.

    Parameters
    ----------
    model, data, stock, rate
        See :func:`residual_table`.
    allow_mc : bool, optional
        Permit the Monte Carlo pricing route.

    Returns
    -------
    tuple[dict[str, float], pd.DataFrame, pd.DataFrame]
        ``(summary, buckets, residuals)``.
    """
    residuals = residual_table(
        model, data, stock, rate, allow_mc=allow_mc, **pricing_kwargs
    )
    return summarize_residuals(residuals), bucket_metrics(residuals), residuals

filter_options

filter_options(df: DataFrame, stock: Stock, rate: Rate, liquidity: LiquidityFilter | None = None) -> pd.DataFrame

Applies a :class:LiquidityFilter to a chain.

Parameters:

Name	Type	Description	Default
df	DataFrame	Raw or annotated chain (see :func:annotate_chain).	required
stock	(Stock, Rate)	Market parameters (forward, implied volatilities).	required
rate	(Stock, Rate)	Market parameters (forward, implied volatilities).	required
liquidity	LiquidityFilter | None	Rules to apply; defaults to LiquidityFilter().	None

Returns:

Type	Description
DataFrame	The surviving options, annotated, with a fresh integer index.

Source code in src/optpricing/calibration/market_filters.py

def filter_options(
    df: pd.DataFrame,
    stock: Stock,
    rate: Rate,
    liquidity: LiquidityFilter | None = None,
) -> pd.DataFrame:
    """
    Applies a :class:`LiquidityFilter` to a chain.

    Parameters
    ----------
    df : pd.DataFrame
        Raw or annotated chain (see :func:`annotate_chain`).
    stock, rate : Stock, Rate
        Market parameters (forward, implied volatilities).
    liquidity : LiquidityFilter | None, optional
        Rules to apply; defaults to ``LiquidityFilter()``.

    Returns
    -------
    pd.DataFrame
        The surviving options, annotated, with a fresh integer index.
    """
    liq = liquidity or LiquidityFilter()
    data = df if "iv_market" in df.columns else annotate_chain(df, stock, rate)
    T = data["maturity"]
    mask = (T >= liq.min_maturity_days / 365.25) & (T <= liq.max_maturity_years)
    mask &= data["moneyness"].between(*liq.moneyness)
    mask &= data["marketPrice"] > 0
    if liq.require_bid and "bid" in data.columns:
        mask &= data["bid"] > 0
    if "relative_spread" in data.columns:
        spread = data["relative_spread"]
        mask &= spread.isna() | (spread <= liq.max_relative_spread)
    if liq.min_open_interest > 0:
        if "openInterest" in data.columns:
            activity = data["openInterest"].fillna(0)
            if "volume" in data.columns:
                activity = activity.where(activity > 0, data["volume"].fillna(0))
            mask &= activity >= liq.min_open_interest
        elif "volume" in data.columns:
            mask &= data["volume"].fillna(0) >= liq.min_open_interest
    if liq.otm_only:
        is_call = data["optionType"] == "call"
        mask &= np.where(is_call, data["moneyness"] >= 1.0, data["moneyness"] < 1.0)
    iv = data["iv_market"]
    mask &= iv.notna() & iv.between(*liq.iv_bounds)
    return data[mask.fillna(False)].reset_index(drop=True)

find_atm_options

find_atm_options(calls: DataFrame, puts: DataFrame, spot: float) -> pd.DataFrame

Finds the closest at-the-money (ATM) call-put pair for each expiry.

Parameters:

Name	Type	Description	Default
calls	DataFrame	A DataFrame of call options.	required
puts	DataFrame	A DataFrame of put options.	required
spot	float	The current spot price of the underlying.	required

Returns:

Type	Description
DataFrame	A DataFrame containing the merged ATM call-put pairs.

Source code in src/optpricing/calibration/fit_market_params.py

def find_atm_options(
    calls: pd.DataFrame,
    puts: pd.DataFrame,
    spot: float,
) -> pd.DataFrame:
    """
    Finds the closest at-the-money (ATM) call-put pair for each expiry.

    Parameters
    ----------
    calls : pd.DataFrame
        A DataFrame of call options.
    puts : pd.DataFrame
        A DataFrame of put options.
    spot : float
        The current spot price of the underlying.

    Returns
    -------
    pd.DataFrame
        A DataFrame containing the merged ATM call-put pairs.
    """
    merged = pd.merge(
        calls, puts, on=["strike", "maturity"], suffixes=("_call", "_put"), how="inner"
    )
    if merged.empty:
        return pd.DataFrame()
    merged["moneyness_dist"] = abs(merged["strike"] - spot)
    atm_indices = merged.groupby("maturity")["moneyness_dist"].idxmin()
    return merged.loc[atm_indices]

fit_jump_params_from_history

fit_jump_params_from_history(log_returns: Series, threshold_stds: float = 3.0) -> dict

Estimates jump parameters and diffusion volatility from historical returns.

This function separates historical log returns into a "diffusion" component (normal daily movements) and a "jump" component (extreme movements) based on a standard deviation threshold. It then calculates the annualized parameters for a jump-diffusion model like Merton's.

Parameters:

Name	Type	Description	Default
log_returns	Series	A pandas Series of daily log returns.	required
threshold_stds	float	The number of standard deviations to use as a threshold for identifying jumps, by default 3.0.	3.0

Returns:

Type	Description
dict	A dictionary containing the estimated parameters: 'sigma', 'lambda', 'mu_j', and 'sigma_j'.

Source code in src/optpricing/calibration/fit_jump_parameters.py

def fit_jump_params_from_history(
    log_returns: pd.Series,
    threshold_stds: float = 3.0,
) -> dict:
    """
    Estimates jump parameters and diffusion volatility from historical returns.

    This function separates historical log returns into a "diffusion" component
    (normal daily movements) and a "jump" component (extreme movements) based
    on a standard deviation threshold. It then calculates the annualized
    parameters for a jump-diffusion model like Merton's.

    Parameters
    ----------
    log_returns : pd.Series
        A pandas Series of daily log returns.
    threshold_stds : float, optional
        The number of standard deviations to use as a threshold for identifying
        jumps, by default 3.0.

    Returns
    -------
    dict
        A dictionary containing the estimated parameters: 'sigma', 'lambda',
        'mu_j', and 'sigma_j'.
    """
    logger.debug("Fitting jump parameters from historical returns...")

    std_dev = log_returns.std()
    jump_threshold = threshold_stds * std_dev

    diffusion_returns = log_returns[abs(log_returns) < jump_threshold]
    jump_returns = log_returns[abs(log_returns) >= jump_threshold]

    # Annualize daily std dev by multiplying by sqrt(252 trading days)
    sigma_est = diffusion_returns.std() * np.sqrt(252)

    if len(jump_returns) > 2:
        lambda_est = len(jump_returns) / len(log_returns) * 252
        mu_j_est = jump_returns.mean()
        sigma_j_est = jump_returns.std()
    else:
        lambda_est, mu_j_est, sigma_j_est = 0.1, 0.0, 0.0
        logger.debug("Not enough jumps detected; using default jump parameters.")

    # Use the key 'sigma' to match the parameter name in the Merton/Kou models.
    fitted_params = {
        "sigma": float(sigma_est),
        "lambda": float(lambda_est),
        "mu_j": float(mu_j_est),
        "sigma_j": float(sigma_j_est),
    }

    formatted = ", ".join(f"{k}: {v:.4f}" for k, v in fitted_params.items())
    logger.debug("Estimated historical params: {%s}", formatted)
    return fitted_params

fit_rate_and_dividend

fit_rate_and_dividend(calls: DataFrame, puts: DataFrame, spot: float, r_fixed: float | None = None, q_fixed: float | None = None) -> tuple[float, float]

Fits the risk-free rate (r) and dividend yield (q) from put-call parity.

With r_fixed given, the dividend yield comes from the implied forward curve (:func:implied_dividend_yield). Otherwise the parity error of the ATM call-put pair of every maturity is minimised over the free parameters (note that parity identifies r - q well but r and q separately only weakly, so supply one of them whenever possible).

Parameters:

Name	Type	Description	Default
calls	DataFrame	A DataFrame of call options with 'strike', 'maturity', 'marketPrice'.	required
puts	DataFrame	A DataFrame of put options with 'strike', 'maturity', 'marketPrice'.	required
spot	float	The current spot price of the underlying.	required
r_fixed	float | None	If provided, the risk-free rate is held fixed at this value. Defaults to None.	None
q_fixed	float | None	If provided, the dividend yield is held fixed at this value. Defaults to None.	None

Returns:

Type	Description
tuple[float, float]	The estimated (or fixed) continuously compounded rate and dividend yield.

Source code in src/optpricing/calibration/fit_market_params.py

def fit_rate_and_dividend(
    calls: pd.DataFrame,
    puts: pd.DataFrame,
    spot: float,
    r_fixed: float | None = None,
    q_fixed: float | None = None,
) -> tuple[float, float]:
    """
    Fits the risk-free rate (r) and dividend yield (q) from put-call parity.

    With ``r_fixed`` given, the dividend yield comes from the implied forward
    curve (:func:`implied_dividend_yield`). Otherwise the parity error of the
    ATM call-put pair of every maturity is minimised over the free parameters
    (note that parity identifies ``r - q`` well but ``r`` and ``q`` separately
    only weakly, so supply one of them whenever possible).

    Parameters
    ----------
    calls : pd.DataFrame
        A DataFrame of call options with 'strike', 'maturity', 'marketPrice'.
    puts : pd.DataFrame
        A DataFrame of put options with 'strike', 'maturity', 'marketPrice'.
    spot : float
        The current spot price of the underlying.
    r_fixed : float | None, optional
        If provided, the risk-free rate is held fixed at this value.
        Defaults to None.
    q_fixed : float | None, optional
        If provided, the dividend yield is held fixed at this value.
        Defaults to None.

    Returns
    -------
    tuple[float, float]
        The estimated (or fixed) continuously compounded rate and dividend yield.
    """
    if r_fixed is not None and q_fixed is not None:
        return float(r_fixed), float(q_fixed)
    if r_fixed is not None:
        return float(r_fixed), implied_dividend_yield(calls, puts, spot, r_fixed)

    atm_pairs = find_atm_options(calls, puts, spot)
    if atm_pairs.empty:
        logger.warning(
            "No call-put pairs found; using r=%.2f, q=%.2f.",
            DEFAULT_RATE,
            DEFAULT_DIVIDEND,
        )
        return DEFAULT_RATE, DEFAULT_DIVIDEND if q_fixed is None else float(q_fixed)

    free_param_indices, initial_guess, bounds = [], [], []
    free_param_indices.append(0)
    initial_guess.append(0.05)
    bounds.append((0.0, 0.15))
    if q_fixed is None:
        free_param_indices.append(1)
        initial_guess.append(0.01)
        bounds.append((0.0, 0.10))

    def parity_error(x: np.ndarray) -> float:
        r = x[0]
        q = q_fixed if 1 not in free_param_indices else x[1]
        parity_rhs = spot * np.exp(-q * atm_pairs["maturity"]) - atm_pairs[
            "strike"
        ] * np.exp(-r * atm_pairs["maturity"])
        parity_lhs = atm_pairs["marketPrice_call"] - atm_pairs["marketPrice_put"]
        return np.sum((parity_lhs - parity_rhs) ** 2)

    solution = minimize(
        fun=parity_error, x0=initial_guess, bounds=bounds, method="L-BFGS-B"
    )
    r_est = float(solution.x[0])
    q_est = float(q_fixed) if q_fixed is not None else float(solution.x[1])
    logger.debug("Fitted market params -> r: %.4f, q: %.4f", r_est, q_est)
    return r_est, q_est

implied_dividend_yield

implied_dividend_yield(calls: DataFrame, puts: DataFrame, spot: float, rate: float, *, min_maturity_days: float = 14.0) -> float

Flat continuous dividend yield implied by the option chain for a given rate.

Aggregates :func:implied_forward_curve with a weighted median (weights T * n_pairs, maturities shorter than min_maturity_days excluded when longer ones exist) and falls back to zero when the result is implausible.

Parameters:

Name	Type	Description	Default
calls	DataFrame	Option quotes (see :func:implied_forward_curve).	required
puts	DataFrame	Option quotes (see :func:implied_forward_curve).	required
spot	float	Spot price.	required
rate	float	Continuously compounded risk-free rate.	required
min_maturity_days	float	Maturities below this are ignored if longer ones exist, by default 14.	14.0

Returns:

Type	Description
float	The implied continuous dividend yield.

Source code in src/optpricing/calibration/fit_market_params.py

def implied_dividend_yield(
    calls: pd.DataFrame,
    puts: pd.DataFrame,
    spot: float,
    rate: float,
    *,
    min_maturity_days: float = 14.0,
) -> float:
    """
    Flat continuous dividend yield implied by the option chain for a given rate.

    Aggregates :func:`implied_forward_curve` with a weighted median (weights
    ``T * n_pairs``, maturities shorter than ``min_maturity_days`` excluded when
    longer ones exist) and falls back to zero when the result is implausible.

    Parameters
    ----------
    calls, puts : pd.DataFrame
        Option quotes (see :func:`implied_forward_curve`).
    spot : float
        Spot price.
    rate : float
        Continuously compounded risk-free rate.
    min_maturity_days : float, optional
        Maturities below this are ignored if longer ones exist, by default 14.

    Returns
    -------
    float
        The implied continuous dividend yield.
    """
    curve = implied_forward_curve(calls, puts, spot, rate)
    if curve.empty:
        return DEFAULT_DIVIDEND
    T = curve["maturity"].to_numpy()
    w = T * curve["n_pairs"].to_numpy()
    long_enough = T >= min_maturity_days / 365.25
    if long_enough.any():
        w = np.where(long_enough, w, 0.0)
    q = _weighted_median(curve["implied_dividend"].to_numpy(), w)
    if not np.isfinite(q) or not (DIVIDEND_BOUNDS[0] <= q <= DIVIDEND_BOUNDS[1]):
        logger.warning("Implied dividend yield %.4f implausible; using 0.", q)
        return DEFAULT_DIVIDEND
    return float(q)

implied_forward_curve

implied_forward_curve(calls: DataFrame, puts: DataFrame, spot: float, rate: float, *, moneyness_window: float = 0.03, min_pairs: int = 3) -> pd.DataFrame

Implied forward (and dividend yield) per maturity from put-call parity.

Each call-put pair gives F_i = K_i + (C_i - P_i) e^{rT}; the maturity's forward is the bid-ask-weighted median over near-the-money pairs. Only pairs close to the money are used because listed US equity and ETF options are American: deep in-the-money puts carry an early-exercise premium that biases parity (it is also why r cannot be inferred from the strike slope of C - P and must be supplied).

Parameters:

Name	Type	Description	Default
calls	DataFrame	Option quotes with strike, maturity, marketPrice and optionally bid/ask columns.	required
puts	DataFrame	Option quotes with strike, maturity, marketPrice and optionally bid/ask columns.	required
spot	float	Spot price used to pick near-the-money pairs and to convert the forward into a dividend yield.	required
rate	float	Continuously compounded risk-free rate.	required
moneyness_window	float	Pairs with |K/S - 1| above this are ignored; the window is widened to the nearest min_pairs strikes when too few qualify. Default 0.03.	0.03
min_pairs	int	Minimum number of pairs per maturity, by default 3.	3

Returns:

Type	Description
DataFrame	One row per maturity: maturity, forward, implied_dividend and n_pairs; empty when no maturity has a usable pair.

Source code in src/optpricing/calibration/fit_market_params.py

def implied_forward_curve(
    calls: pd.DataFrame,
    puts: pd.DataFrame,
    spot: float,
    rate: float,
    *,
    moneyness_window: float = 0.03,
    min_pairs: int = 3,
) -> pd.DataFrame:
    """
    Implied forward (and dividend yield) per maturity from put-call parity.

    Each call-put pair gives ``F_i = K_i + (C_i - P_i) e^{rT}``; the maturity's
    forward is the bid-ask-weighted median over near-the-money pairs. Only
    pairs close to the money are used because listed US equity and ETF options
    are American: deep in-the-money puts carry an early-exercise premium that
    biases parity (it is also why ``r`` cannot be inferred from the strike
    slope of ``C - P`` and must be supplied).

    Parameters
    ----------
    calls, puts : pd.DataFrame
        Option quotes with ``strike``, ``maturity``, ``marketPrice`` and
        optionally ``bid``/``ask`` columns.
    spot : float
        Spot price used to pick near-the-money pairs and to convert the
        forward into a dividend yield.
    rate : float
        Continuously compounded risk-free rate.
    moneyness_window : float, optional
        Pairs with ``|K/S - 1|`` above this are ignored; the window is widened
        to the nearest ``min_pairs`` strikes when too few qualify. Default 0.03.
    min_pairs : int, optional
        Minimum number of pairs per maturity, by default 3.

    Returns
    -------
    pd.DataFrame
        One row per maturity: ``maturity``, ``forward``, ``implied_dividend``
        and ``n_pairs``; empty when no maturity has a usable pair.
    """
    cols = ["strike", "maturity", "marketPrice"]
    has_quotes = all(c in calls.columns and c in puts.columns for c in ("bid", "ask"))
    quote_cols = ["bid", "ask"] if has_quotes else []
    merged = pd.merge(
        calls[cols + quote_cols],
        puts[cols + quote_cols],
        on=["strike", "maturity"],
        suffixes=("_call", "_put"),
        how="inner",
    ).dropna(subset=["marketPrice_call", "marketPrice_put"])
    if merged.empty:
        return pd.DataFrame()

    rows = []
    for T, grp in merged.groupby("maturity"):
        if T <= 0:
            continue
        grp = grp.drop_duplicates("strike")
        dist = (grp["strike"] / spot - 1.0).abs()
        near = grp[dist <= moneyness_window]
        if len(near) < min_pairs:
            near = grp.loc[dist.nsmallest(min_pairs).index]
        if near.empty:
            continue
        K = near["strike"].to_numpy(dtype=float)
        cmp_ = (near["marketPrice_call"] - near["marketPrice_put"]).to_numpy(float)
        forwards = K + cmp_ * math.exp(rate * T)
        if has_quotes:
            spread = (
                near["ask_call"] - near["bid_call"] + near["ask_put"] - near["bid_put"]
            ).to_numpy(dtype=float)
            w = 1.0 / np.maximum(spread, 1e-3)
        else:
            w = np.ones_like(K)
        forward = _weighted_median(forwards, w)
        if not forward > 0:
            continue
        rows.append(
            {
                "maturity": float(T),
                "forward": float(forward),
                "implied_dividend": float(rate - math.log(forward / spot) / T),
                "n_pairs": int(len(near)),
            }
        )
    return pd.DataFrame(rows)

price_options_vectorized

price_options_vectorized(options_df: DataFrame, stock: Stock, model: BaseModel, rate: Rate, *, upper_bound: float = 1000.0, method: str = 'cos', allow_mc: bool = False, mc_paths: int = 20000, mc_steps: int = 50, mc_seed: int = 0, **kwargs: Any) -> np.ndarray

Prices a chain of European options with the fastest exact method available.

Parameters:

Name	Type	Description	Default
options_df	DataFrame	Must contain strike, maturity and optionType. The index is reset internally to guarantee safe positional writes.	required
stock	Stock	Underlying description.	required
model	BaseModel	A model with a vectorized closed form, a characteristic function or (with allow_mc) a simulation method.	required
rate	Rate	Continuous zero-curve.	required
upper_bound	float	Integration truncation limit of the method="gauss" route.	1000
method	('cos', 'gauss')	Characteristic-function pricer (see :class:ChainLayout).	"cos"
allow_mc	bool	Fall back to Monte Carlo for simulation-only models.	False
mc_paths	int	Monte Carlo settings for that fallback (fixed seed: common random numbers across calls).	20000
mc_steps	int	Monte Carlo settings for that fallback (fixed seed: common random numbers across calls).	20000
mc_seed	int	Monte Carlo settings for that fallback (fixed seed: common random numbers across calls).	20000
**kwargs	Any	Extra pricing arguments forwarded to the characteristic function or simulation kernel (e.g. v0).	{}

Returns:

Type	Description
ndarray	Model prices - aligned with the row order of options_df.

Raises:

Type	Description
NotImplementedError	If the model has neither a vectorized closed form nor a characteristic function and allow_mc is False.

Source code in src/optpricing/calibration/vectorized_pricer.py

def price_options_vectorized(
    options_df: pd.DataFrame,
    stock: Stock,
    model: BaseModel,
    rate: Rate,
    *,
    upper_bound: float = 1000.0,
    method: str = "cos",
    allow_mc: bool = False,
    mc_paths: int = 20_000,
    mc_steps: int = 50,
    mc_seed: int = 0,
    **kwargs: Any,
) -> np.ndarray:
    """
    Prices a chain of European options with the fastest exact method available.

    Parameters
    ----------
    options_df : pd.DataFrame
        Must contain `strike`, `maturity` and `optionType`.
        The index is reset internally to guarantee safe positional writes.
    stock : Stock
        Underlying description.
    model : BaseModel
        A model with a vectorized closed form, a characteristic function or
        (with ``allow_mc``) a simulation method.
    rate : Rate
        Continuous zero-curve.
    upper_bound : float, default 1000
        Integration truncation limit of the ``method="gauss"`` route.
    method : {"cos", "gauss"}, default "cos"
        Characteristic-function pricer (see :class:`ChainLayout`).
    allow_mc : bool, default False
        Fall back to Monte Carlo for simulation-only models.
    mc_paths, mc_steps, mc_seed : int
        Monte Carlo settings for that fallback (fixed seed: common random
        numbers across calls).
    **kwargs : Any
        Extra pricing arguments forwarded to the characteristic function or
        simulation kernel (e.g. ``v0``).

    Returns
    -------
    np.ndarray
        Model prices - aligned with the row order of options_df.

    Raises
    ------
    NotImplementedError
        If the model has neither a vectorized closed form nor a characteristic
        function and ``allow_mc`` is False.
    """
    layout = ChainLayout(options_df, upper_bound=upper_bound, method=method)
    return price_chain(
        layout,
        stock,
        model,
        rate,
        allow_mc=allow_mc,
        mc_paths=mc_paths,
        mc_steps=mc_steps,
        mc_seed=mc_seed,
        **kwargs,
    )

residual_table

residual_table(model: BaseModel, data: DataFrame, stock: Stock, rate: Rate, *, allow_mc: bool = False, iv_solver: BSMIVSolver | None = None, **pricing_kwargs: Any) -> pd.DataFrame

Prices every option in data with model and returns the residuals.

Parameters:

Name	Type	Description	Default
model	BaseModel	The (calibrated) model.	required
data	DataFrame	Options with strike, maturity, optionType, marketPrice; annotated columns (iv_market, moneyness) are computed if absent.	required
stock	(Stock, Rate)	Market parameters.	required
rate	(Stock, Rate)	Market parameters.	required
allow_mc	bool	Permit the Monte Carlo pricing route, by default False.	False
iv_solver	BSMIVSolver | None	Solver for the model implied volatilities.	None
**pricing_kwargs	Any	Forwarded to :func:price_options_vectorized.	{}

Returns:

Type	Description
DataFrame	data plus model_price, iv_model, price_error (model - market), rel_price_error and iv_error.

Source code in src/optpricing/calibration/metrics.py

def residual_table(
    model: BaseModel,
    data: pd.DataFrame,
    stock: Stock,
    rate: Rate,
    *,
    allow_mc: bool = False,
    iv_solver: BSMIVSolver | None = None,
    **pricing_kwargs: Any,
) -> pd.DataFrame:
    """
    Prices every option in ``data`` with ``model`` and returns the residuals.

    Parameters
    ----------
    model : BaseModel
        The (calibrated) model.
    data : pd.DataFrame
        Options with ``strike``, ``maturity``, ``optionType``, ``marketPrice``;
        annotated columns (``iv_market``, ``moneyness``) are computed if absent.
    stock, rate : Stock, Rate
        Market parameters.
    allow_mc : bool, optional
        Permit the Monte Carlo pricing route, by default False.
    iv_solver : BSMIVSolver | None, optional
        Solver for the model implied volatilities.
    **pricing_kwargs : Any
        Forwarded to :func:`price_options_vectorized`.

    Returns
    -------
    pd.DataFrame
        ``data`` plus ``model_price``, ``iv_model``, ``price_error``
        (model - market), ``rel_price_error`` and ``iv_error``.
    """
    table = data if "iv_market" in data.columns else annotate_chain(data, stock, rate)
    table = table.reset_index(drop=True).copy()
    if table.empty:
        for col in ("model_price", "iv_model", "price_error", "rel_price_error"):
            table[col] = np.nan
        table["iv_error"] = np.nan
        return table
    prices = price_options_vectorized(
        table, stock, model, rate, allow_mc=allow_mc, **pricing_kwargs
    )
    solver = iv_solver or BSMIVSolver()
    iv_model = solver.solve(
        prices, table, stock, rate, initial_guess=table["iv_market"].to_numpy()
    )
    market = table["marketPrice"].to_numpy(dtype=float)
    table["model_price"] = prices
    table["iv_model"] = iv_model
    table["price_error"] = prices - market
    with np.errstate(divide="ignore", invalid="ignore"):
        table["rel_price_error"] = (prices - market) / market
    table["iv_error"] = iv_model - table["iv_market"].to_numpy()
    return table

select_fastest_technique

select_fastest_technique(model: BaseModel)

Selects the fastest available pricing technique for a given model.

The function checks the capabilities of the model in a specific order of preference, which generally corresponds to computational speed.

The order of preference is: 1. Closed-Form 2. COS (Fang & Oosterlee 2008), for characteristic-function models: as accurate as FFT and 0.1-0.5 ms per option against 0.4-55 ms 3. Monte Carlo

Parameters:

Name	Type	Description	Default
model	BaseModel	The financial model for which to select a technique.	required

Returns:

Type	Description
BaseTechnique	An instance of the fastest suitable pricing technique.

Raises:

Type	Description
TypeError	If no suitable pricing technique can be found for the model.

Source code in src/optpricing/calibration/technique_selector.py

def select_fastest_technique(model: BaseModel):
    """
    Selects the fastest available pricing technique for a given model.

    The function checks the capabilities of the model in a specific order of
    preference, which generally corresponds to computational speed.

    The order of preference is:
    1. Closed-Form
    2. COS (Fang & Oosterlee 2008), for characteristic-function models: as
       accurate as FFT and 0.1-0.5 ms per option against 0.4-55 ms
    3. Monte Carlo

    Parameters
    ----------
    model : BaseModel
        The financial model for which to select a technique.

    Returns
    -------
    BaseTechnique
        An instance of the fastest suitable pricing technique.

    Raises
    ------
    TypeError
        If no suitable pricing technique can be found for the model.
    """
    if model.has_closed_form:
        return ClosedFormTechnique()
    if model.supports_cf:
        return COSTechnique()
    if model.supports_sde:
        return MonteCarloTechnique(n_paths=5000, n_steps=50, antithetic=True)
    raise TypeError(f"No suitable pricing technique found for model '{model.name}'")

summarize_residuals

summarize_residuals(table: DataFrame) -> dict[str, float]

Scalar metrics from a :func:residual_table.

Returns:

Type	Description
dict[str, float]	n_options, rmse_price, mae_price, rmse_rel_price, rmse_iv (vol points as a fraction), mae_iv, n_iv_failed (options whose model price has no Black-Scholes volatility).

Source code in src/optpricing/calibration/metrics.py

def summarize_residuals(table: pd.DataFrame) -> dict[str, float]:
    """
    Scalar metrics from a :func:`residual_table`.

    Returns
    -------
    dict[str, float]
        ``n_options``, ``rmse_price``, ``mae_price``, ``rmse_rel_price``,
        ``rmse_iv`` (vol points as a fraction), ``mae_iv``, ``n_iv_failed``
        (options whose model price has no Black-Scholes volatility).
    """
    if table.empty:
        return {
            "n_options": 0,
            "rmse_price": float("nan"),
            "mae_price": float("nan"),
            "rmse_rel_price": float("nan"),
            "rmse_iv": float("nan"),
            "mae_iv": float("nan"),
            "n_iv_failed": 0,
        }
    iv_err = table["iv_error"].to_numpy(dtype=float)
    return {
        "n_options": int(len(table)),
        "rmse_price": _rmse(table["price_error"]),
        "mae_price": float(np.nanmean(np.abs(table["price_error"]))),
        "rmse_rel_price": _rmse(table["rel_price_error"]),
        "rmse_iv": _rmse(iv_err),
        "mae_iv": float(np.nanmean(np.abs(iv_err)))
        if np.isfinite(iv_err).any()
        else float("nan"),
        "n_iv_failed": int(np.sum(~np.isfinite(iv_err))),
    }

supports_vectorized_pricing

supports_vectorized_pricing(model: BaseModel, *, allow_mc: bool = False) -> bool

True if :func:price_options_vectorized can price a chain for model.

Parameters:

Name	Type	Description	Default
model	BaseModel	The model to check.	required
allow_mc	bool	Count simulation-only models as supported, by default False.	False

Source code in src/optpricing/calibration/vectorized_pricer.py

def supports_vectorized_pricing(model: BaseModel, *, allow_mc: bool = False) -> bool:
    """
    True if :func:`price_options_vectorized` can price a chain for ``model``.

    Parameters
    ----------
    model : BaseModel
        The model to check.
    allow_mc : bool, optional
        Count simulation-only models as supported, by default False.
    """
    if model.has_closed_form and getattr(model, "vectorized_closed_form", False):
        return True
    if model.supports_cf:
        return True
    return allow_mc and (
        model.supports_sde
        or getattr(model, "is_pure_levy", False)
        or getattr(model, "has_exact_sampler", False)
    )