用代码理解导数

它变化得有多快？

函数 f 在 x 处的导数是它的瞬时变化率：曲线在该点的斜率。在金融中，几乎所有敏感度都是导数。期权的 Delta 是标的每变动一个单位其价格变动多少，Gamma 是 Delta 本身变化得多快，债券久期是价格对收益率的导数。微积分为简单函数给出精确公式（x**2 的导数是 2*x），但真实的定价函数往往太复杂，或只以代码形式存在。这时我们用数值方法来估计导数。

# A function is just a Python function
def f(x):
    return x * x

# Slope between two nearby points approaches the true slope (2x = 6 at x = 3)
x = 3.0
for h in (1.0, 0.1, 0.001):
    slope = (f(x + h) - f(x)) / h
    print(h, round(slope, 4))
# 1.0 7.0
# 0.1 6.1
# 0.001 6.001

有限差分：前向与中心

前向差分向前看一小步 h：(f(x + h) - f(x)) / h，误差与 h 成正比地减小。中心差分向两侧各看一步：(f(x + h) - f(x - h)) / (2 * h)，两侧的误差相互抵消，误差按 h 的平方减小，对二次函数是精确的。二阶导数（Gamma）结合三个点：(f(x + h) - 2 * f(x) + f(x - h)) / h**2。把函数本身作为参数传入（f），这是前面函数章节学过的技巧。

import math

def forward_diff(f, x, h=1e-6):
    return (f(x + h) - f(x)) / h

def central_diff(f, x, h=1e-5):
    return (f(x + h) - f(x - h)) / (2 * h)

def second_derivative(f, x, h=1e-4):
    return (f(x + h) - 2 * f(x) + f(x - h)) / h**2

print(round(central_diff(math.sin, 0.0), 6))            # 1.0  (cos 0)
print(round(central_diff(lambda x: x**3, 2.0), 4))      # 12.0 (3 x^2)
print(round(second_derivative(lambda x: x**3, 2.0), 3)) # 12.0 (6 x)

# A quadratic: central is exact even for a big step, forward is not
print(central_diff(lambda x: x * x, 3.0, h=0.5))  # 6.0
print(forward_diff(lambda x: x * x, 3.0, h=0.5))  # 6.5

选择步长

h 越小，公式的截断误差越小，但两个几乎相等的数相减会损失有效数字：这种舍入误差随 h 变小而增大，二阶导数还要除以 h**2，增长更快。h 太小（如 1e-12）可能得到毫无意义的结果。折中的好选择是：一阶导数用中心差分时约 1e-5，二阶导数约 1e-4。在真实的定价代码里，“bump and reprice”（扰动后重新定价）说的就是这个：把现货价格小幅上移和下移，各重新运行一次定价器。

import math

def central_diff(f, x, h):
    return (f(x + h) - f(x - h)) / (2 * h)

# True derivative of exp at 0 is 1; watch the error as h shrinks
for h in (1e-1, 1e-3, 1e-5, 1e-8, 1e-12):
    error = abs(central_diff(math.exp, 0.0, h) - 1.0)
    print(f"h={h:g} error={error:.1e}")
# h=0.1 error=1.7e-03
# h=0.001 error=1.7e-07
# h=1e-05 error=1.2e-11
# h=1e-08 error=5.3e-10     (rounding error takes over)
# h=1e-12 error=3.3e-05     (much worse: too small)

导数的应用：牛顿法

导数告诉你朝哪个方向走才能到达根，也就是 f(x) = 0 的点。牛顿法从一个猜测值 x 出发，反复用 x - f(x) / f'(x) 替换它，沿着切线一路走向零点。它在根附近收敛得非常快，因此常被用来由期权价格反推隐含波动率。当 abs(f(x)) 小于容差时停止，超过最大迭代次数就放弃；如果导数为零或循环不收敛，应当报错，而不是返回一个错误的数。

def central_diff(f, x, h=1e-5):
    return (f(x + h) - f(x - h)) / (2 * h)

def newton_root(f, x0, tol=1e-10, max_iter=50):
    x = x0
    for _ in range(max_iter):
        if abs(f(x)) < tol:
            return x
        slope = central_diff(f, x)
        if slope == 0:
            raise ValueError("flat slope: Newton cannot continue")
        x = x - f(x) / slope
    raise ValueError("did not converge")

# The square root of 2 is the positive root of x^2 - 2
root = newton_root(lambda x: x * x - 2, 1.0)
print(round(root, 8))           # 1.41421356

# x^2 + 1 never reaches zero
try:
    newton_root(lambda x: x * x + 1, 1.0)
except ValueError:
    print("no root found")      # no root found

• 前向 (f(x+h) - f(x)) / h；中心 (f(x+h) - f(x-h)) / (2h)；二阶 (f(x+h) - 2f(x) + f(x-h)) / h**2。
• 不大不小：一阶约 1e-5，二阶约 1e-4。
• 牛顿法：重复 x = x - f(x) / f'(x) 直到 abs(f(x)) 很小；不收敛则报错。

下一节：用代码求积分，把许多细条加起来，求面积、总量和期望值。

实验：评分内容

• 用 forward_diff 向前估计斜率25 分
• 用 central_diff 从两侧估计斜率25 分
• 用 second_derivative 估计曲率（Gamma）25 分
• 用 newton_root 找到函数为零的位置25 分