模型校准

让模型与市场一致

定价模型有参数，正确的取值不在教科书里，而在今天的市场里。校准就是选取参数，使模型价格尽可能接近观察到的价格，通常是当天的期权报价。做法有三部分。一个衡量偏差的目标函数，几乎总是差值的平方和，可以加权以便流动性好的报价占更大比重。一个移动参数以减小它的搜索。以及对结果是否合理的检验：逐日稳定、在合理范围内，并且在没参与拟合的报价上表现良好。

# The objective: how far the model's prices are from the market's, as a (weighted) sum of squares
def sum_squared_error(model, market, weights=None):
    if weights is None:
        weights = [1.0] * len(model)
    return sum(w * (m - y) ** 2 for w, m, y in zip(weights, model, market))

market = [4.20, 8.60, 14.10]
model = [4.00, 8.80, 14.20]
print(round(sum_squared_error(model, market), 4))                    # 0.09   0.2^2 + 0.2^2 + 0.1^2
print(round(sum_squared_error(model, market, [1, 4, 1]), 4))         # 0.21   the middle quote counts four times as much

两种搜索工具：求根与求最小值

有两项任务不断出现。第一是隐含波动率：给定市场价格，找出使模型给出该价格的波动率，也就是 model(sigma) - price = 0 的根。因为价格随波动率单调上升，二分法可行：从差值变号的区间出发，对半切开，保留符号仍变化的那一半，然后重复；每一步精度翻倍，而且不会失败。第二是最小值：对只有一个谷的函数，黄金分割搜索每步把区间缩小到原来的 0.618，并重用两个函数值之一，所以每步只需一次求值，也不需要导数。库（scipy.optimize）更快并能处理很多参数，但这两个工具能展示思路。

# Implied volatility by bisection: the root of (model price - market price) in sigma
import math

def norm_cdf(x):
    return 0.5 * (1 + math.erf(x / math.sqrt(2)))

def bs_call(s, k, t, r, sigma):
    d1 = (math.log(s / k) + (r + sigma ** 2 / 2) * t) / (sigma * math.sqrt(t))
    return s * norm_cdf(d1) - k * math.exp(-r * t) * norm_cdf(d1 - sigma * math.sqrt(t))

def bisect_root(f, lo, hi, tol=1e-10):
    f_lo = f(lo)
    while hi - lo > tol:
        mid = (lo + hi) / 2
        f_mid = f(mid)
        if f_lo * f_mid <= 0:
            hi = mid
        else:
            lo, f_lo = mid, f_mid
    return (lo + hi) / 2

market_price = bs_call(100.0, 105.0, 1.0, 0.02, 0.27)           # a quote generated with a volatility of 27%
implied = bisect_root(lambda sigma: bs_call(100.0, 105.0, 1.0, 0.02, sigma) - market_price, 0.01, 2.0)
print(round(implied, 6))                                         # 0.27   recovered

# Golden-section search: one function evaluation per step, the interval shrinks by 0.618 each time
import math

def golden_section_minimize(f, lo, hi, tol=1e-8):
    ratio = (math.sqrt(5) - 1) / 2
    a, b = lo, hi
    c, d = b - ratio * (b - a), a + ratio * (b - a)
    fc, fd = f(c), f(d)
    while b - a > tol:
        if fc < fd:
            b, d, fd = d, c, fc
            c = b - ratio * (b - a)
            fc = f(c)
        else:
            a, c, fc = c, d, fd
            d = a + ratio * (b - a)
            fd = f(d)
    return (a + b) / 2

evaluations = []
def valley(x):
    evaluations.append(x)
    return (x - 1.3) ** 2 + 0.5

print(round(golden_section_minimize(valley, 0.0, 5.0), 5), len(evaluations) < 60)       # 1.3 True   about 40 evaluations for 1e-8

校准一个波动率，并且负责任地做

把这些部件组合起来，把单个波动率校准到多个期权价格上，就是对 sigma 最小化 price_fn(strike, sigma) 与市场价格之间的平方误差之和，price_fn 作为参数传入，所以同一段代码适用于任何模型。对于有多个参数的模型（Heston 有五个，SABR 固定 beta 后有三个），搜索是多维的，可能陷入局部极小：从多个初值出发，给参数设边界，并优先使用相对误差或隐含波动率误差而不是原始价格误差，以免忽略便宜的虚值期权。然后检查拟合出的参数逐日是否稳定，以及模型对没参与拟合的期权定价是否准确。完美却不稳定的拟合就是过拟合。

# Calibrate one volatility to three option prices: minimise the sum of squared errors over sigma
import math

def norm_cdf(x):
    return 0.5 * (1 + math.erf(x / math.sqrt(2)))

def bs_call(s, k, t, r, sigma):
    d1 = (math.log(s / k) + (r + sigma ** 2 / 2) * t) / (sigma * math.sqrt(t))
    return s * norm_cdf(d1) - k * math.exp(-r * t) * norm_cdf(d1 - sigma * math.sqrt(t))

def golden_section_minimize(f, lo, hi, tol=1e-8):
    ratio = (math.sqrt(5) - 1) / 2
    a, b = lo, hi
    c, d = b - ratio * (b - a), a + ratio * (b - a)
    fc, fd = f(c), f(d)
    while b - a > tol:
        if fc < fd:
            b, d, fd = d, c, fc
            c = b - ratio * (b - a)
            fc = f(c)
        else:
            a, c, fc = c, d, fd
            d = a + ratio * (b - a)
            fd = f(d)
    return (a + b) / 2

def calibrate_sigma(strikes, market_prices, price_fn, lo, hi):
    def error(sigma):
        return sum((price_fn(k, sigma) - p) ** 2 for k, p in zip(strikes, market_prices))
    return golden_section_minimize(error, lo, hi)

strikes = [90.0, 100.0, 110.0]
quotes = [bs_call(100.0, k, 1.0, 0.02, 0.23) for k in strikes]                  # market quotes made with sigma = 23%
sigma = calibrate_sigma(strikes, quotes, lambda k, s: bs_call(100.0, k, 1.0, 0.02, s), 0.01, 1.5)
print(round(sigma, 4))                                                           # 0.23

• 校准选取使价格误差（加权）平方和最小的模型参数；要检查逐日稳定性以及在未参与拟合的报价上的精度。
• 二分法通过把变号区间对半切分来求根（隐含波动率）；黄金分割用每步一次求值找出单谷函数的最小值。
• 把定价函数作为参数传入，给参数设界，多参数模型从多个初值出发，并拒绝荒谬的拟合。

下一节：用于预测的神经网络，从数据中学习预测，而不是自己写公式。

实验：评分内容

• 用 sum_squared_error 衡量偏差25 分
• 用 bisect_root 求解 f(x) = 025 分
• 用 golden_section_minimize 找最低点25 分
• 用 calibrate_sigma 拟合市场25 分