返回 Expert 笔记
Expert Day 67

Week 10 复习 — options_lib v0.1 整合

复习 Day 61-66:概率/随机过程/BS/Greeks/IV/SABR 的关系

2026-07-07
Phase 2 - 量化数学与衍生品定价 (Day 61-74)
量化期权整合单元测试options_lib

日期: 2026-07-07 方向: 量化 / 衍生品定价 阶段: Phase 2 - 量化数学与衍生品定价 (Day 61-74) 标签: #量化 #期权 #整合 #单元测试 #options_lib


今日目标

类型内容
学习复习 Day 61-66:概率/随机过程/BS/Greeks/IV/SABR 的关系
实操整合代码为 options_lib v0.1,加单元测试,统一 API
产出完整 options_lib/ 包 + 测试 + Cookbook

一、知识地图

Day 61 概率统计基础 ─────┐
                         ├─→ Day 62 随机过程 (BM, Itô, GBM)
                         │       ↓
                         │   Day 63 BS PDE 与封闭解
                         │       ↓
                         ├─→ Day 64 Greeks (一阶/二阶)
                         │       ↓
                         ├─→ Day 65 IV 反解 + Smile
                         │       ↓
                         └─→ Day 66 SABR 随机波动率

核心关联

  • Day 61 的 CLT 是 Day 62 BM 的"无穷小极限"
  • Day 62 的 GBM SDE 在 Day 63 应用 Itô 引理推 BS PDE
  • Day 64 的 Greeks 是 BS 公式对参数的偏导
  • Day 65 IV 反解暴露 BS 的不足(smile)
  • Day 66 SABR 是修正 BS 缺陷的动态模型

二、整合架构 — options_lib v0.1

2.1 包结构

options_lib/
├── __init__.py
├── core.py          # BS 公式 + Greeks
├── implied.py       # IV 反解 (Newton, Brent)
├── smile.py         # SVI, SABR 模型
├── greeks.py        # 二阶 Greeks
├── stochastic.py    # GBM/SDE 模拟 (来自 Day 62)
├── stats.py         # 矩估计、CLT 工具 (来自 Day 61)
└── tests/
    ├── test_core.py
    ├── test_implied.py
    └── test_smile.py

2.2 统一 API 设计

from options_lib import EuropeanOption, ImpliedVolSolver, SABR

# 1. 基本定价
opt = EuropeanOption(S=65000, K=70000, r=0.05, sigma=0.65, tau=30/365, q=0, kind="call")
price = opt.price()
greeks = opt.greeks()  # 返回 dict

# 2. IV 反解
iv = ImpliedVolSolver.solve(market_price=1654.0, S=65000, K=70000, r=0.05, tau=30/365, kind="call")

# 3. SABR 校准
sabr = SABR(beta=1.0)
sabr.calibrate(F=65000, T=30/365, strikes=[...], market_ivs=[...])
print(sabr.alpha, sabr.rho, sabr.nu)
sigma_at_k = sabr.iv(K=80000)  # 在 K=80000 的 SABR-implied IV

三、完整代码:options_lib/

3.1 core.py

"""options_lib/core.py - BS 公式与一阶 Greeks"""

import numpy as np
from scipy.stats import norm
from dataclasses import dataclass, field
from typing import Literal


@dataclass
class EuropeanOption:
    S: float                    # 现货价
    K: float                    # 行权价
    r: float                    # 无风险利率 (年化)
    sigma: float                # 波动率 (年化)
    tau: float                  # 距到期 (年)
    q: float = 0.0              # 股息/staking yield
    kind: Literal["call", "put"] = "call"

    def __post_init__(self):
        assert self.S > 0 and self.K > 0
        assert self.tau >= 0 and self.sigma > 0

    @property
    def d1(self):
        return (np.log(self.S/self.K) + (self.r - self.q + 0.5*self.sigma**2)*self.tau) \
               / (self.sigma * np.sqrt(self.tau))

    @property
    def d2(self):
        return self.d1 - self.sigma * np.sqrt(self.tau)

    def price(self):
        if self.tau <= 0:
            return max(self.S - self.K, 0) if self.kind == "call" else max(self.K - self.S, 0)
        d1, d2 = self.d1, self.d2
        S, K, r, q, tau = self.S, self.K, self.r, self.q, self.tau
        if self.kind == "call":
            return S*np.exp(-q*tau)*norm.cdf(d1) - K*np.exp(-r*tau)*norm.cdf(d2)
        else:
            return K*np.exp(-r*tau)*norm.cdf(-d2) - S*np.exp(-q*tau)*norm.cdf(-d1)

    def delta(self):
        if self.kind == "call":
            return np.exp(-self.q*self.tau) * norm.cdf(self.d1)
        else:
            return np.exp(-self.q*self.tau) * (norm.cdf(self.d1) - 1)

    def gamma(self):
        return np.exp(-self.q*self.tau) * norm.pdf(self.d1) / (self.S*self.sigma*np.sqrt(self.tau))

    def vega(self, per_pct=True):
        v = self.S * np.exp(-self.q*self.tau) * norm.pdf(self.d1) * np.sqrt(self.tau)
        return v / 100 if per_pct else v

    def theta(self, per_day=True):
        d1, d2 = self.d1, self.d2
        S, K, r, q, sigma, tau = self.S, self.K, self.r, self.q, self.sigma, self.tau
        common = -S*np.exp(-q*tau)*norm.pdf(d1)*sigma / (2*np.sqrt(tau))
        if self.kind == "call":
            t = common - r*K*np.exp(-r*tau)*norm.cdf(d2) + q*S*np.exp(-q*tau)*norm.cdf(d1)
        else:
            t = common + r*K*np.exp(-r*tau)*norm.cdf(-d2) - q*S*np.exp(-q*tau)*norm.cdf(-d1)
        return t / 365 if per_day else t

    def rho(self, per_pct=True):
        if self.kind == "call":
            r_ = self.K * self.tau * np.exp(-self.r*self.tau) * norm.cdf(self.d2)
        else:
            r_ = -self.K * self.tau * np.exp(-self.r*self.tau) * norm.cdf(-self.d2)
        return r_ / 100 if per_pct else r_

    def greeks(self):
        return {
            "price": self.price(),
            "delta": self.delta(),
            "gamma": self.gamma(),
            "vega":  self.vega(),
            "theta": self.theta(),
            "rho":   self.rho(),
        }

3.2 implied.py

"""options_lib/implied.py - IV 反解"""

import numpy as np
from scipy.optimize import brentq
from .core import EuropeanOption


class ImpliedVolSolver:

    @staticmethod
    def solve_newton(market_price, S, K, r, tau, q=0.0, kind="call",
                     sigma0=None, tol=1e-7, max_iter=100):
        if sigma0 is None:
            sigma0 = max(0.1, np.sqrt(2*abs(np.log(S/K) + (r-q)*tau)/tau))
        sigma = sigma0
        for _ in range(max_iter):
            opt = EuropeanOption(S, K, r, sigma, tau, q, kind)
            diff = opt.price() - market_price
            if abs(diff) < tol:
                return sigma
            v = opt.vega(per_pct=False)
            if v < 1e-12:
                break
            sigma -= diff / v
            sigma = max(1e-4, min(5.0, sigma))
        return sigma

    @staticmethod
    def solve_brent(market_price, S, K, r, tau, q=0.0, kind="call",
                    sigma_lo=1e-4, sigma_hi=5.0):
        def f(s):
            return EuropeanOption(S, K, r, s, tau, q, kind).price() - market_price
        try:
            return brentq(f, sigma_lo, sigma_hi, xtol=1e-7)
        except ValueError:
            return np.nan

    @staticmethod
    def solve(market_price, S, K, r, tau, q=0.0, kind="call", method="newton"):
        if method == "newton":
            return ImpliedVolSolver.solve_newton(market_price, S, K, r, tau, q, kind)
        elif method == "brent":
            return ImpliedVolSolver.solve_brent(market_price, S, K, r, tau, q, kind)
        else:
            raise ValueError(f"Unknown method: {method}")

3.3 smile.py

"""options_lib/smile.py - SVI / SABR 模型"""

import numpy as np
from scipy.optimize import minimize
from dataclasses import dataclass


@dataclass
class SABR:
    beta: float = 1.0
    alpha: float = None
    rho: float = None
    nu: float = None
    F: float = None
    T: float = None

    def iv(self, K):
        return _hagan_iv(self.F, K, self.T, self.alpha, self.beta, self.rho, self.nu)

    def calibrate(self, F, T, strikes, market_ivs, x0=None):
        self.F, self.T = F, T
        strikes = np.asarray(strikes); market_ivs = np.asarray(market_ivs)
        if x0 is None:
            atm_idx = np.argmin(np.abs(strikes - F))
            x0 = [market_ivs[atm_idx] * F**(1-self.beta), 0.0, 0.4]
        def loss(p):
            a, r, n = p
            if a <= 0 or abs(r) >= 1 or n < 0:
                return 1e10
            ivs = np.array([_hagan_iv(F, K, T, a, self.beta, r, n) for K in strikes])
            if np.any(np.isnan(ivs)):
                return 1e10
            return np.sum((ivs - market_ivs)**2)
        bounds = [(1e-4, 5), (-0.999, 0.999), (1e-4, 5)]
        result = minimize(loss, x0, bounds=bounds, method="L-BFGS-B")
        self.alpha, self.rho, self.nu = result.x
        return self


def _hagan_iv(F, K, T, alpha, beta, rho, nu):
    if abs(F - K) < 1e-10:
        FK_b = F**(1-beta)
        A = ((1-beta)**2/24)*alpha**2/FK_b**2 \
            + 0.25*rho*beta*nu*alpha/FK_b \
            + (2-3*rho**2)/24*nu**2
        return alpha / FK_b * (1 + A*T)
    log_FK = np.log(F/K)
    FK_g = (F*K)**((1-beta)/2)
    z = (nu/alpha) * FK_g * log_FK
    if abs(z) < 1e-10:
        zx = 1.0
    else:
        sq = np.sqrt(1 - 2*rho*z + z*z)
        x_ = np.log((sq + z - rho)/(1-rho))
        zx = z / x_
    denom = 1 + ((1-beta)**2/24)*log_FK**2 + ((1-beta)**4/1920)*log_FK**4
    main = alpha / (FK_g * denom)
    A = ((1-beta)**2/24)*alpha**2/FK_g**2 \
        + 0.25*rho*beta*nu*alpha/FK_g \
        + (2-3*rho**2)/24*nu**2
    return main * zx * (1 + A*T)


@dataclass
class SVI:
    a: float = None
    b: float = None
    rho: float = None
    m: float = None
    s: float = None
    F: float = None
    T: float = None

    def total_variance(self, k):
        return self.a + self.b*(self.rho*(k-self.m) + np.sqrt((k-self.m)**2 + self.s**2))

    def iv(self, K):
        k = np.log(K/self.F)
        w = self.total_variance(k)
        return np.sqrt(np.maximum(w, 1e-8) / self.T)

    def calibrate(self, F, T, strikes, market_ivs):
        self.F, self.T = F, T
        k = np.log(np.asarray(strikes) / F)
        w_market = (np.asarray(market_ivs)**2) * T
        def loss(p):
            a, b, r, m, s = p
            w = a + b*(r*(k-m) + np.sqrt((k-m)**2 + s**2))
            return np.sum((w - w_market)**2)
        bounds = [(-1, 1), (1e-6, 5), (-0.999, 0.999), (-2, 2), (1e-6, 2)]
        x0 = [np.median(w_market), 0.1, -0.3, 0.0, 0.1]
        result = minimize(loss, x0, bounds=bounds, method="L-BFGS-B")
        self.a, self.b, self.rho, self.m, self.s = result.x
        return self

3.4 tests/test_core.py

"""单元测试"""
import numpy as np
import pytest
from options_lib.core import EuropeanOption
from options_lib.implied import ImpliedVolSolver
from options_lib.smile import SABR


def test_put_call_parity():
    """C - P = S - K e^(-rT)"""
    opt_c = EuropeanOption(S=100, K=100, r=0.05, sigma=0.2, tau=1.0, kind="call")
    opt_p = EuropeanOption(S=100, K=100, r=0.05, sigma=0.2, tau=1.0, kind="put")
    diff = opt_c.price() - opt_p.price()
    parity = 100 - 100 * np.exp(-0.05 * 1.0)
    assert abs(diff - parity) < 1e-10


def test_atm_delta_half():
    """ATM Call Delta 接近 0.5 (短期, 0 利率)"""
    opt = EuropeanOption(S=100, K=100, r=0.0, sigma=0.2, tau=0.001, kind="call")
    assert abs(opt.delta() - 0.5) < 0.01


def test_call_vs_intrinsic_value():
    """Call 价 >= max(S - K, 0)"""
    for K in [80, 100, 120]:
        opt = EuropeanOption(S=100, K=K, r=0.05, sigma=0.3, tau=0.5, kind="call")
        assert opt.price() >= max(100 - K, 0) - 1e-6


def test_iv_round_trip():
    """反解 IV: BS 价 -> IV -> BS 价应一致"""
    opt = EuropeanOption(S=100, K=110, r=0.05, sigma=0.25, tau=0.5, kind="call")
    price = opt.price()
    iv = ImpliedVolSolver.solve(price, 100, 110, 0.05, 0.5, kind="call")
    assert abs(iv - 0.25) < 1e-5


def test_sabr_atm_consistency():
    """SABR 在 K=F 处应给出 ATM IV"""
    sabr = SABR(beta=1.0)
    sabr.alpha, sabr.rho, sabr.nu = 0.3, -0.2, 0.5
    sabr.F, sabr.T = 100, 1.0
    iv_at_F = sabr.iv(100)
    assert 0.25 < iv_at_F < 0.35


def test_gamma_positivity():
    """Gamma 始终 >= 0"""
    for K in [80, 100, 120]:
        opt = EuropeanOption(S=100, K=K, r=0.05, sigma=0.3, tau=0.5, kind="call")
        assert opt.gamma() >= 0


if __name__ == "__main__":
    import sys
    pytest.main([__file__, "-v"])

3.5 Cookbook 示例

"""使用示例 cookbook.py"""
from options_lib import EuropeanOption, ImpliedVolSolver, SABR
import numpy as np

# ===== 1. 基础定价 =====
opt = EuropeanOption(S=65000, K=70000, r=0.05, sigma=0.65, tau=30/365, kind="call")
print(f"Price: ${opt.price():.2f}")
print(f"Greeks: {opt.greeks()}")

# ===== 2. IV 反解 =====
mkt_price = 1654.36
iv = ImpliedVolSolver.solve(mkt_price, S=65000, K=70000, r=0.05, tau=30/365, kind="call")
print(f"Implied Vol: {iv*100:.2f}%")

# ===== 3. SABR 校准 =====
F = 65000
T = 30/365
strikes = [50000, 55000, 60000, 65000, 70000, 75000, 80000]
ivs = [0.68, 0.61, 0.56, 0.55, 0.58, 0.64, 0.72]

sabr = SABR(beta=1.0)
sabr.calibrate(F, T, strikes, ivs)
print(f"SABR: alpha={sabr.alpha:.4f}, rho={sabr.rho:+.4f}, nu={sabr.nu:.4f}")

# 在 K=72500 处 interpolate IV
iv_72500 = sabr.iv(72500)
print(f"SABR IV at K=72500: {iv_72500*100:.2f}%")

# 用 SABR-IV 给 BTC-72500-C 定价
opt_72500 = EuropeanOption(S=F, K=72500, r=0.05, sigma=iv_72500, tau=T, kind="call")
print(f"BTC-72500-C price (SABR-IV): ${opt_72500.price():.2f}")

四、单元测试运行

$ pytest options_lib/tests/ -v
============================= test session starts ==============================
collected 6 items

test_core.py::test_put_call_parity PASSED                                [ 16%]
test_core.py::test_atm_delta_half PASSED                                 [ 33%]
test_core.py::test_call_vs_intrinsic_value PASSED                        [ 50%]
test_core.py::test_iv_round_trip PASSED                                  [ 66%]
test_core.py::test_sabr_atm_consistency PASSED                           [ 83%]
test_core.py::test_gamma_positivity PASSED                               [100%]

============================== 6 passed in 0.42s ===============================

五、加密市场特化整合

options_lib 中专门加 crypto-specific 模块:

"""options_lib/crypto.py - 加密市场专用"""

class DeribitInverseOption(EuropeanOption):
    """
    Deribit BTC/ETH inverse 期权:
    premium 以 BTC 报价, 而非 USD
    """
    def price_btc(self):
        """以 BTC 计价"""
        return self.price() / self.S

    def delta_inverse(self):
        """Inverse Delta (premium 单位 = base currency)"""
        delta_usd = self.delta()
        c_usd = self.price()
        return delta_usd / self.S - c_usd / self.S**2


class PerpetualOption:
    """Power Perp (Squeeth-style) 期权"""
    def __init__(self, S, alpha, n_power=2):
        self.S = S
        self.alpha = alpha  # 瞬时方差
        self.n = n_power

    def funding_rate(self):
        """Power Perp funding = (n^2 - n)/2 * sigma^2 (per Squeeth paper)"""
        return self.n*(self.n - 1) / 2 * self.alpha

六、常见陷阱(整合版)

  1. 类型一致:用 dataclass 让参数显式且可校验。

  2. 数值稳定:所有 sigma/tau 必须 > 0;IV 反解加 bounds 防跑飞。

  3. API 一致性:所有定价对象返回相同 dict 字段(price/delta/gamma/vega/theta/rho)。

  4. 测试覆盖:parity、ATM 极限、round-trip、Greeks 正负号必须测。

  5. 加密扩展:inverse delta、power perp 单独继承类。


七、关键速查

模块内容
core.pyEuropeanOption, BS price + 一阶 Greeks
implied.pyImpliedVolSolver.solve(method='newton'/'brent')
smile.pySABR(beta=1).calibrate(...), SVI().calibrate(...)
crypto.pyDeribitInverseOption, PerpetualOption
tests/pytest 单元测试
命令用途
pytest options_lib/tests/跑单元测试
mypy options_lib/类型检查

八、面试题

Q1: 设计期权定价库的 API,你会怎么权衡 OO vs 函数式?

主接口用 OO(EuropeanOption(...).price())让对象封装参数;底层数学函数用纯函数(便于向量化、测试、JIT 优化)。混合最优。

Q2: 怎么测试期权定价库的正确性?

(1) 已知封闭解对照(教科书数值);(2) Put-Call Parity;(3) Greeks 一致性(finite difference vs 解析);(4) MC 收敛验证;(5) IV round-trip。

Q3: 生产环境跑 SABR 校准 1000 个 IV 表面(不同 underlying × 不同 expiry),怎么优化?

(1) 校准并行化(每个 surface 独立);(2) 用 JAX/numba JIT 加速 Hagan 公式;(3) warm-start:今天用昨天的参数作初值;(4) 缓存 N(d) 计算(向量化所有 strikes)。

Q4: 库里 BS 公式遇到 $\sigma = 0$ 或 $\tau = 0$ 怎么处理?

单独 branch:$\tau \to 0$ 返回 intrinsic value;$\sigma = 0$ 返回 forward 与 strike 比较的 deterministic payoff(即 $e^{-rT} \max(F - K, 0)$)。

Q5: 如果用户报告 BS 反解 IV 在某些深 OTM 期权"卡住",怎么调试?

(1) 检查 market price < intrinsic value(套利机会,IV 不存在);(2) Vega 接近 0,Newton 失败 → 用 Brent;(3) market price > upper bound($S e^{-qT}$ for Call),同样无解;(4) 加 sanity check 抛出明确错误。


九、明日预告

Day 68: 美式期权 — 二叉树定价、有限差分、最优行权边界。明天我们离开欧式期权封闭解,进入需要数值方法的世界。Cox-Ross-Rubinstein (CRR) 二叉树是最经典的离散方法,并能优雅处理早行权特征。