6.17. 求解器与随机数

当所研究的函数由 Python 代码而非样本缓冲区定义时,另一类工具才是合适的选择:函数的根在哪里、它的最小值是多少、它在给定区间上的积分是多少?scipy.integratescipy.optimize 子模块涵盖了这类工作。每个算法都会回调到用户提供的 Python 函数中,因此每次迭代的开销比缓冲区归约要高;其便利之处在于无需自己编写求解器。

scipy.special 子模块涵盖了统计专用函数(误差函数、伽马函数),在计算概率分布的累积分布函数(CDF,即样本不超过给定值的概率)或概率密度函数(PDF,即在给定值处的相对似然)时会用到它们。numpy.random 涵盖了用于抖动、仿真和合成测试数据的伪随机生成器。

6.17.1. 可调用对象的数值积分

当被积函数是一个 Python 函数而非样本缓冲区时,scipy.integrate 提供了四种求积算法:

  • quad() —— 自适应 Gauss-Kronrod 法。对于光滑被积函数,这是合适的默认选择。返回 (value, error)

  • romberg() —— 经典的 Romberg / Newton-Cotes 法。返回单个浮点数。上游已弃用;为兼容性而保留。

  • simpson() —— 自适应 Simpson 法则。返回单个浮点数。

  • tanhsinh() —— 双指数求积法。当被积函数在端点存在奇点或具有无限积分限时使用。返回 (value, error)

用双指数法则(tanhsinh)计算的高斯积分:

from math import exp
from math import pi
from math import sqrt
from ulab import numpy as np
from ulab import scipy as sp

f = lambda x: exp(-x * x)
value, err = sp.integrate.tanhsinh(f, -np.inf, np.inf)
print("approx:", value, "   exact:", sqrt(pi))

输出:

approx: 1.7724538...   exact: 1.7724538...

6.17.2. 求根与最小化

scipy.optimize 涵盖了三种经典的单变量求解器。每次迭代都会回调到用户提供的 Python 函数中,因此相比纯 Python 求解器的加速是有限的(大约 2 倍);其便利之处在于无需自己编写求解器。

  • bisect() —— 通过对区间二分法在 [a, b] 上查找 f 的根。f(a)f(b) 必须符号相反:

    def f(x):
        return x * x - 1
    
    sp.optimize.bisect(f, 0, 4)        # ~1.0
    
  • newton() —— 使用割线法 / Newton-Raphson 迭代查找根:

    def f(x):
        return x * x * x - 2.0
    
    sp.optimize.newton(f, 3., tol=0.001, rtol=0.01)
    # ~1.260
    
  • fmin() —— 使用下降单纯形(Nelder-Mead)法查找局部最小值:

    def f(x):
        return (x - 1) ** 2 - 1
    
    sp.optimize.fmin(f, 3.0)           # ~1.0
    

单变量的适用范围足以满足大多数摄像头侧的优化——传感器的标定常数、使对比度测量最大化的增益、直方图双峰性最显著处的阈值。对于多变量问题,正确的做法通常是将问题重新表述为一个小型线性代数求解,而不是诉诸通用的非线性优化器。

6.17.3. 专用函数

scipy.special 提供了少数几个行为类似通用函数的统计与概率函数——它们接受标量、可迭代对象或 ndarray,并返回一个浮点 ndarray:

x = np.linspace(0, 4, num=8)

sp.special.erf(x)         # error function
sp.special.erfc(x)        # complementary error function
sp.special.gamma(x + 1)   # gamma function
sp.special.gammaln(x + 1) # log-gamma function

误差函数及其补函数出现在高斯分布的 CDF 中——它们是在已测得的 z 分数与概率之间进行转换,或计算正态分布尾部积分时的首选方法。伽马函数和对数伽马函数出现在 beta / 卡方 / 学生 t 计算中;gammaln 是大参数下的数值稳定形式,此时 gamma 本身会溢出。

6.17.4. 随机数

numpy.random 提供了一个 Generator 类,从常见分布中抽取样本。该生成器是有状态的:每次调用都会推进其内部状态,因此连续的调用会返回相互独立的样本:

from ulab import numpy as np

rng = np.random.Generator(seed=42)

rng.random(size=5)             # 5 uniform [0.0, 1.0) samples
rng.uniform(low=-1.0, high=1.0, size=10)
rng.normal(loc=0.0, scale=1.0, size=(2, 4))

输出 dtype 始终为 floatsize= 接受一个整数(一维输出)或一个元组(n 维输出);省略时返回单个 Python 浮点数。

该生成器适用于仿真、抖动、合成测试数据,以及任何不要求加密强度的应用。它适用于密钥或令牌;这类用途请通过 os 使用系统随机源。

6.17.5. 构建时可用性

每个子模块是否实际存在取决于摄像头的构建方式。scipy.optimizescipy.special 并非在每台摄像头上都启用;调用摄像头未包含的函数会引发 AttributeErrordir(sp)dir(sp.optimize)dir(np.random) 等可报告目标摄像头上有哪些可用项。