尧图网站建设 
PyMC 协方差函数Covariance Functions完全指南pymc.gp.cov 模块深度解析【免费下载链接】pymcBayesian Modeling and Probabilistic Programming in Python项目地址: https://gitcode.com/GitHub_Trending/py/pymc高斯过程Gaussian Process, GP的核心在于协方差函数又称核函数k(x, x)它决定了函数先验的平滑性、周期性与各向异性等一切性质。本指南以 PyMC 官方 API 文档 docs/source/api/gp/cov.rst 为骨架结合 pymc/gp/cov.py 源码与 tests/gp/test_cov.py 测试用例系统讲解 PyMC 中pymc.gp.cov模块的全部内置协方差函数、input_dim与active_dims的参数语义、核的代数组合加、乘、标量缩放、幂运算以及如何在gp.Latent、gp.Marginal等 GP 实现中落地使用。读完本文你将能够针对不同数据特征选择并组合出合适的协方差函数构建可解释、可组合的贝叶斯 GP 模型。协方差函数在 PyMC GP 中的角色在贝叶斯建模中有时我们关心的未知量不是标量或定长向量而是一个连续函数。高斯过程为函数f(x)提供了先验分布f(x) ~ GP(m(x), k(x, x))其中m(x)是均值函数k(x, x)是协方差函数。函数值被建模为多元正态分布的一次抽样其协方差结构完全由k(x, x)决定。PyMC 的 GP 模块正是利用多元正态分布的边缘化与条件化性质来完成推断边缘分布与预测条件分布的。PyMC 的 GP 具备清晰语法 高度可组合两大特点参见官方指南 docs/source/guides/Gaussian_Processes.rst模块内置了大量预定义的协方差函数、均值函数和多种 GP 实现更重要的是GP 在 PyMC 中被当作可嵌入更大层次模型中的分布而非只能独立使用的回归器。pymc.gp.cov正是全部协方差函数的实现模块。根据 API 文档该模块通过automodule:: pymc.gp.cov自动生成文档公开的类包括Constant、WhiteNoise、ExpQuad、RatQuad、Exponential、Matern52、Matern32、Linear、Polynomial、Cosine、Periodic、WarpedInput、Gibbs、Coregion、ScaledCov、Kron。此外源码__all__pymc/gp/cov.py中还包含Matern12与WrappedPeriodic。实例化与求值先参数化后喂数据使用过 GPy 或 GPflow 的用户会对 PyMC 的语法感到熟悉协方差函数在实例化时只完成参数化并不接触输入数据真正求值发生在之后以cov_func(X, Xs)方式调用时。这一延迟绑定设计正是为了支持核的组合构造——组合出的新核在被调用前无需关心具体输入。所有协方差函数都继承自基类BaseCovariancepymc/gp/cov.py其统一的调用接口为cov_func(X, XsNone, diagFalse)X训练输入形状(n, input_dim)Xs可选的预测输入若为None则等价于Xs Xdiag为True时只返回协方差矩阵的对角元用于降低计算量例如仅需方差时。调用内部会分发到两个抽象方法BaseCovariancefull(X, Xs)计算完整的n x n或n x m协方差矩阵diag(X)仅计算对角元素。以最简单的核为例import pymc as pm import numpy as np X np.linspace(0, 1, 10)[:, None] # 10 个一维输入点 cov pm.gp.cov.ExpQuad(input_dim1, ls0.1) K cov(X).eval() # 10x10 协方差矩阵 Kd cov(X, diagTrue).eval() # 只取对角input_dim 与 active_dims告诉核操作哪些列除Constant和WhiteNoise外绝大多数协方差函数继承自Covariance基类pymc/gp/cov.py其构造函数统一接收两个关键参数参数含义默认值input_dim输入矩阵X的总列数总输入维度必填active_dims该核实际操作的列索引列表全部列np.arange(input_dim)之所以必须显式声明input_dim是因为协方差函数在构造时尚未看到任何输入数据。active_dims允许同一个核只作用于输入的部分维度从而可以构建不同核作用于不同维度、再组合的模型。官方指南给出了经典示例docs/source/guides/Gaussian_Processes.rst对包含三个预测变量的矩阵构造一个只作用于第二、三列的指数二次核并为每个维度配置独立的长度尺度ls [2, 5] # 第二、三列各自的长度尺度 cov_func pm.gp.cov.ExpQuad(input_dim3, lsls, active_dims[1, 2])源码中_slice方法pymc/gp/cov.py会按active_dims对X和Xs做列切片并在此处给出一个实用警告若X的实际列数与input_dim不一致会提示只有input_dim列被用于计算协方差提醒用户确认参数意图。需要特别说明的参数校验Covariance.initactive_dims中的任何值都不能超过input_dim否则抛出ValueError。内置协方差函数逐一详解以下按照 API 文档列出的顺序结合源码中的数学定义与参数说明逐个展开。常数与噪声核Constantpymc/gp/cov.pyk(x, x) c恒为常数c的协方差函数。它不继承Covariance无input_dim/active_dims常作为偏置项参与核组合。求值时通过_alloc生成全c矩阵。WhiteNoisepymc/gp/cov.pyk(x, x) σ²·I白噪声协方差函数sigma为噪声标准差。其对角为σ²当给定预测输入Xs时full返回n x m的全零矩阵——这正确反映了不同输入点之间白噪声不相关的性质。常被加在核上建模观测噪声cov_func ExpQuad(...) WhiteNoise(sigma)。平稳核Stationary所有平稳核继承Stationary基类pymc/gp/cov.py其共性是核值只依赖输入点之间的距离因此diag恒为 1.0。Stationary统一处理长度尺度参数构造函数签名如下Stationary(input_dim, lsNone, ls_invNone, active_dimsNone)ls长度尺度。input_dim 1时可以是标量列表/数组每个维度一个也支持 PyMC 随机变量input_dim 1时为标量。ls_inv逆长度尺度即1 / ls。ls与ls_inv必须且只能提供一个否则抛出ValueErrorpymc/gp/cov.py。在源码中ls_inv会被转换为ls 1.0 / ls_inv。距离计算方面square_dist先按1/ls缩放输入再计算平方距离并用pt.clip(sqd, 0.0, np.inf)消除数值误差导致的负值euclidean_dist在平方距离上开方并加上1e-12数值稳定项pymc/gp/cov.py。ExpQuad指数二次核又称平方指数 / RBF 核pymc/gp/cov.pyk(x, x) exp( -(x - x)² / (2ℓ²) )最常用的光滑核处处无穷可微。full_from_distance实现为pt.exp(-0.5 * r2)。它还实现了功率谱密度见后文谱密度一节。RatQuad有理二次核pymc/gp/cov.pyk(x, x) ( 1 (x - x)² / (2αℓ²) )^(-α)额外参数alphaα控制形状它可以看作无穷多个不同长度尺度的平方指数核的尺度混合α 越大越接近ExpQuad。源码注释给出了其作为平方指数核关于精度参数 λ ~ Gamma(α, αℓ²) 的混合表示。Exponential指数核pymc/gp/cov.pyk(x, x) exp( -||x - x|| / (2ℓ) )与 ExpQuad 不同它基于欧氏距离而非平方距离对应的是奥恩斯坦-乌伦贝克OU过程样本路径连续但不可导。Matern52 / Matern32 / Matern12Matérn 核族Matern52ν 5/2pymc/gp/cov.pyk(x, x) (1 √5·r/ℓ 5r²/(3ℓ²)) · exp(-√5·r/ℓ)Matern32ν 3/2pymc/gp/cov.pyk(x, x) (1 √3·r/ℓ) · exp(-√3·r/ℓ)Matern12ν 1/2pymc/gp/cov.pyk(x, x) exp(-r/ℓ)其中r为经长度尺度缩放后的欧氏距离。Matérn 族通过 ν 参数刻画函数可微性ν 越大函数越光滑Matern52二阶可微、Matern32一阶可微、Matern12即指数核的变体不可微。实践中Matern52是平衡光滑性与数值稳定性的常用选择。Cosine余弦核pymc/gp/cov.pyk(x, x) cos( 2π·||x - x|| / ℓ² )强周期振荡核适用于周期性但非平滑衰减的信号。Periodic周期核pymc/gp/cov.pyk(x, x) exp( -sin²(π|x - x|/T) / (2ℓ²) )需要额外参数period周期T。源码特别给出了一个易踩坑的说明PyMC 的系数约定与常见定义不同指数上是0.5而非常见的2因此当你想复现标准定义时初始化时需将长度尺度除以 2pymc/gp/cov.py。此外Periodic还实现了power_spectral_density_approx用于 HSGP 低秩近似的系数计算基于第一类修正贝塞尔函数I_j。非平稳核Linear线性核pymc/gp/cov.pyk(x, x) (x - c)(x - c)c为偏移中心点。该核产生的函数是输入空间的线性函数常用于趋势建模。Polynomial多项式核pymc/gp/cov.pyk(x, x) [ (x - c)(x - c) offset ]^d在Linear基础上增加次数d与常数项offset实现方式为对线性核的结果取(linear offset)^d。输入变换类核WarpedInputpymc/gp/cov.pyk(x, x) k_base(w(x), w(x))用任意 PyTensor 函数warp_func对输入做非线性扭曲再交给底层核cov_func计算。参数args用于向warp_func传递额外的标量或 PyMC 变量内部由handle_args包装见 pymc/gp/cov.py。典型应用将一维输入映射到高维特征空间从而获得非平稳行为。Gibbspymc/gp/cov.pyk(x, x) sqrt( 2ℓ(x)ℓ(x) / (ℓ²(x) ℓ²(x)) ) · exp( -(x - x)² / (ℓ²(x) ℓ²(x)) )使用随输入变化的长度尺度函数lengthscale_func构造非平稳核。源码明确标注仅在 1 维下测试过若active_dims长度大于 1 或input_dim ! 1会直接抛出NotImplementedErrorpymc/gp/cov.py。args同样用于传递额外参数。ScaledCovpymc/gp/cov.pyk(x, x) φ(x) · k_base(x, x) · φ(x)用非负的缩放函数φ(x)scaling_funcPyTensor 可调用对象对基础核cov_func做点级缩放实现对函数幅值的非平稳调制。源码中diag实现为cov_diag * φ(X)²full实现为outer(φ(X), φ(Xs)) * k_base(X, Xs)。多输出核CoregionCoregion协区域化核pymc/gp/cov.py用于内禀/线性协区域化模型ICM/LCM其协方差矩阵为B W·Wᵀ diag(κ)W形状(num_outputs, rank)的低秩矩阵决定各输出之间的相关性kappa形状(num_outputs,)的向量使各输出可独立变化B形状(num_outputs, num_outputs)的完整矩阵。约束条件(W, kappa)与B二者必须恰好提供一个并且该核要求恰好激活一个维度active_dims长度必须为 1否则抛出ValueErrorpymc/gp/cov.py。调用时输入应为整数索引输出编号源码将其cast为int32后通过B[index, index2]查表取值。结构核KronKronKronecker 积核pymc/gp/cov.py与普通乘法所有核共享同一份输入不同Kron 核先把输入按各因子的input_dim切分到各自子空间再对各子空间分别求核矩阵最后做 Kronecker 乘积Kron([cov1, cov2, ...])其整体input_dim等于各因子input_dim之和因子必须是协方差函数或其组合数组不支持。它通常配合gp.MarginalKron与gp.LatentKron使用用于输入为笛卡尔积结构的数据如网格数据可显著加速矩阵运算。测试中也使用pymc.math.kronecker来验证其正确性tests/gp/test_cov.py。核的代数组合可组合性是 PyMC GP 的灵魂协方差函数在 PyMC 中严格遵循核函数的代数规则——这一点被设计为语言级特性BaseCovariance 中的运算符重载用户可以像拼接乐高一样构造复杂核两个核相加仍是核Addpymc/gp/cov.pycov_func pm.gp.cov.ExpQuad(input_dim1, ls1.0) pm.gp.cov.Periodic(input_dim1, period0.5)两个核相乘仍是核Prodpymc/gp/cov.pycov_func pm.gp.cov.ExpQuad(input_dim1, ls1.0) * pm.gp.cov.Periodic(input_dim1, period0.5)核与标量的乘/加仍是核标量会自动被包装为Constantcov_func eta**2 * pm.gp.cov.Matern32(input_dim1, ls1.0) # 幅值缩放 cov_func pm.gp.cov.ExpQuad(input_dim1, ls1.0) 1.0 # 加常数偏置幂运算Exponentiatedpymc/gp/cov.pycov_func ** p要求底数必须是继承自Covariance的核Constant/WhiteNoise会报TypeError且指数必须是标量。底层实现上Combination基类pymc/gp/cov.py做了两件事自动推导元信息组合核的input_dim取各因子之交集所有因子必须一致否则报错active_dims取各因子激活维度的并集排序后。扁平化因子列表Add套Add、Prod套Prod会被自动展平避免嵌套过深。_merge_factors_cov还负责统一处理diagTrue时对数组/张量因子取对角np.diag/pt.diag的细节保证对角计算的正确性——这一点被测试用例反复验证tests/gp/test_cov.py。从测试可以看出这些运算的实际数值行为tests/gp/test_cov.pyExpQuad(1, 0.1) 1的非对角元为1.53940正是ExpQuad核值0.53940加上常数1ExpQuad(1, 0.1) WhiteNoise(sigma1)的对角元为2核对角 1 加噪声方差 1、非对角元保持0.53940。功率谱密度从频域理解平稳核对于平稳核PyMC 还实现了power_spectral_density(omega)Stationary为 HSGPHilbert Space Gaussian Process等谱方法提供理论基础ExpQuad的谱密度为高斯形式pymc/gp/cov.pyRatQuad的谱密度涉及第二类修正贝塞尔函数K_ν(z)并处理了z 0处的奇异性pymc/gp/cov.pyMatern32/Matern52的谱密度为有理函数形式pymc/gp/cov.py。组合核的谱密度遵循简单规则Add的谱密度等于各因子谱密度之和Prod仅当因子中协方差函数不多于一个时才能计算两个协方差函数相乘的谱没有解析形式会抛出NotImplementedError见 pymc/gp/cov.py。此外进行谱密度计算时要求求和项中所有核的active_dims完全一致pymc/gp/cov.py非平稳核则会收到明确的ValueError提示。在 GP 模型中使用协方差函数协方差函数本身只是零件真正发挥威力的是与 GP 实现的配合。PyMC 提供gp.Latent、gp.Marginal等多种实现统一的使用模式为先用均值函数与协方差函数实例化 GP 对象再调用prior、marginal_likelihood或conditional方法构造 PyMC 随机变量。Latent GP 的完整流程docs/source/guides/Gaussian_Processes.rstimport pymc as pm # 1. 定义核与均值 cov_func pm.gp.cov.ExpQuad(input_dim1, ls1.0) mean_func pm.gp.mean.Zero() # 2. 实例化 GP gp pm.gp.Latent(mean_func, cov_func) with pm.Model() as model: # 3. 在观测点 X 上建立函数先验 f f gp.prior(f, X) # 4. 预测点的条件分布 f* f_star gp.conditional(f_star, X_star)prior的第一个参数是随机变量名第二个是函数输入X通常来自数据也可以是 PyMC 随机变量——若输入是张量/随机变量则必须显式给出shape。gp.Marginal类没有prior方法改用marginal_likelihood需要额外提供观测数据与噪声等参数。可加 GPAdditive GP是 PyMC 的另一大特性docs/source/guides/Gaussian_Processes.rstGP 对象之间可以直接相加从而把复杂函数分解为多个独立成分并分别推断with pm.Model() as model: gp1 pm.gp.Marginal(mean_func1, cov_func1) # 长程趋势 gp2 pm.gp.Marginal(mean_func2, cov_func2) # 周期成分 gp gp1 gp2 # f1 f2 f gp.marginal_likelihood(f, X, y, noise) idata pm.sample(1000)注意相加的两个 GP 对象类型必须一致gp.Marginal不能与gp.Latent相加。若想得到某个成分如gp2的条件分布需要在conditional中以given{X: X, y: y, noise: noise, gp: gp}字典形式显式传入缓存的参数docs/source/guides/Gaussian_Processes.rstwith model: f2_star gp2.conditional(f2_star, X_star, given{X: X, y: y, noise: noise, gp: gp}) f_star gp.conditional(f_star, X_star) # 整体条件分布无需 given整体gp的条件分布无需given因为参数已在marginal_likelihood调用时被缓存而gp1、gp2自身从未调用过marginal_likelihood必须显式提供。测试驱动的工程细节与注意事项仓库中的 tests/gp/test_cov.py共 948 行是对pymc.gp.cov行为的最权威验证值得注意的工程细节包括对角一致性所有测试都校验cov(X, diagTrue)与np.diag(cov(X))严格一致如 tests/gp/test_cov.py这是BaseCovariance接口设计正确性的基本保证。标量/矩阵的左右结合a cov、cov a、M cov、cov M均被支持通过__radd__/__rmul__与__array_wrap__但形状非法的组合如三维数组会抛出ValueError: cannot combine...tests/gp/test_cov.py。组合核的元信息合并不同input_dim的核相加会报错active_dims取并集这保证组合核始终能正确切片输入。实战中的几条核心建议长度尺度选择ls可设为标量各向同性或按维度设置的数组自动相关性确定ARD当维度较多时优先用ls_inv配合正约束先验如HalfNormal、Gamma进行推断。噪声与核分离观测噪声建议显式使用WhiteNoise或gp.Marginal的noise参数保持核本身干净便于解释。周期数据优先选用Periodic注意其 0.5 缩放约定若底层核需要更强的可微性可用WrappedPeriodic将任意Stationary核卷成周期核pymc/gp/cov.py。大规模数据网格化/可分解结构数据用Kron配合gp.MarginalKron/gp.LatentKron否则可结合Periodic.power_spectral_density_approx与 HSGP 近似见 pymc/gp/hsgp_approx.py降低计算负担。小结pymc.gp.cov为贝叶斯高斯过程建模提供了完整、可组合、经过充分测试的协方差函数工具箱从经典的ExpQuad、Matérn 族、周期核到非平稳的WarpedInput、Gibbs、ScaledCov再到多输出Coregion与结构化的Kron全部遵循统一的BaseCovariance调用契约并支持加、乘、标量缩放与幂运算的代数组合。配合input_dim/active_dims的维度控制与ls/ls_inv的长度尺度参数化你可以在 pymc/gp/gp.py 提供的Latent、Marginal、MarginalKron等实现中自由搭建从简单回归到可加成分模型的各种 GP 模型。【免费下载链接】pymcBayesian Modeling and Probabilistic Programming in Python项目地址: https://gitcode.com/GitHub_Trending/py/pymc创作声明：本文部分内容由AI辅助生成（AIGC），仅供参考 PyMC 协方差函数（Covariance Functions）完全指南：pymc.gp.cov 模块深度解析
PyMC 协方差函数Covariance Functions完全指南pymc.gp.cov 模块深度解析【免费下载链接】pymcBayesian Modeling and Probabilistic Programming in Python项目地址: https://gitcode.com/GitHub_Trending/py/pymc高斯过程Gaussian Process, GP的核心在于协方差函数又称核函数k(x, x)它决定了函数先验的平滑性、周期性与各向异性等一切性质。本指南以 PyMC 官方 API 文档 docs/source/api/gp/cov.rst 为骨架结合 pymc/gp/cov.py 源码与 tests/gp/test_cov.py 测试用例系统讲解 PyMC 中pymc.gp.cov模块的全部内置协方差函数、input_dim与active_dims的参数语义、核的代数组合加、乘、标量缩放、幂运算以及如何在gp.Latent、gp.Marginal等 GP 实现中落地使用。读完本文你将能够针对不同数据特征选择并组合出合适的协方差函数构建可解释、可组合的贝叶斯 GP 模型。协方差函数在 PyMC GP 中的角色在贝叶斯建模中有时我们关心的未知量不是标量或定长向量而是一个连续函数。高斯过程为函数f(x)提供了先验分布f(x) ~ GP(m(x), k(x, x))其中m(x)是均值函数k(x, x)是协方差函数。函数值被建模为多元正态分布的一次抽样其协方差结构完全由k(x, x)决定。PyMC 的 GP 模块正是利用多元正态分布的边缘化与条件化性质来完成推断边缘分布与预测条件分布的。PyMC 的 GP 具备清晰语法 高度可组合两大特点参见官方指南 docs/source/guides/Gaussian_Processes.rst模块内置了大量预定义的协方差函数、均值函数和多种 GP 实现更重要的是GP 在 PyMC 中被当作可嵌入更大层次模型中的分布而非只能独立使用的回归器。pymc.gp.cov正是全部协方差函数的实现模块。根据 API 文档该模块通过automodule:: pymc.gp.cov自动生成文档公开的类包括Constant、WhiteNoise、ExpQuad、RatQuad、Exponential、Matern52、Matern32、Linear、Polynomial、Cosine、Periodic、WarpedInput、Gibbs、Coregion、ScaledCov、Kron。此外源码__all__pymc/gp/cov.py中还包含Matern12与WrappedPeriodic。实例化与求值先参数化后喂数据使用过 GPy 或 GPflow 的用户会对 PyMC 的语法感到熟悉协方差函数在实例化时只完成参数化并不接触输入数据真正求值发生在之后以cov_func(X, Xs)方式调用时。这一延迟绑定设计正是为了支持核的组合构造——组合出的新核在被调用前无需关心具体输入。所有协方差函数都继承自基类BaseCovariancepymc/gp/cov.py其统一的调用接口为cov_func(X, XsNone, diagFalse)X训练输入形状(n, input_dim)Xs可选的预测输入若为None则等价于Xs Xdiag为True时只返回协方差矩阵的对角元用于降低计算量例如仅需方差时。调用内部会分发到两个抽象方法BaseCovariancefull(X, Xs)计算完整的n x n或n x m协方差矩阵diag(X)仅计算对角元素。以最简单的核为例import pymc as pm import numpy as np X np.linspace(0, 1, 10)[:, None] # 10 个一维输入点 cov pm.gp.cov.ExpQuad(input_dim1, ls0.1) K cov(X).eval() # 10x10 协方差矩阵 Kd cov(X, diagTrue).eval() # 只取对角input_dim 与 active_dims告诉核操作哪些列除Constant和WhiteNoise外绝大多数协方差函数继承自Covariance基类pymc/gp/cov.py其构造函数统一接收两个关键参数参数含义默认值input_dim输入矩阵X的总列数总输入维度必填active_dims该核实际操作的列索引列表全部列np.arange(input_dim)之所以必须显式声明input_dim是因为协方差函数在构造时尚未看到任何输入数据。active_dims允许同一个核只作用于输入的部分维度从而可以构建不同核作用于不同维度、再组合的模型。官方指南给出了经典示例docs/source/guides/Gaussian_Processes.rst对包含三个预测变量的矩阵构造一个只作用于第二、三列的指数二次核并为每个维度配置独立的长度尺度ls [2, 5] # 第二、三列各自的长度尺度 cov_func pm.gp.cov.ExpQuad(input_dim3, lsls, active_dims[1, 2])源码中_slice方法pymc/gp/cov.py会按active_dims对X和Xs做列切片并在此处给出一个实用警告若X的实际列数与input_dim不一致会提示只有input_dim列被用于计算协方差提醒用户确认参数意图。需要特别说明的参数校验Covariance.initactive_dims中的任何值都不能超过input_dim否则抛出ValueError。内置协方差函数逐一详解以下按照 API 文档列出的顺序结合源码中的数学定义与参数说明逐个展开。常数与噪声核Constantpymc/gp/cov.pyk(x, x) c恒为常数c的协方差函数。它不继承Covariance无input_dim/active_dims常作为偏置项参与核组合。求值时通过_alloc生成全c矩阵。WhiteNoisepymc/gp/cov.pyk(x, x) σ²·I白噪声协方差函数sigma为噪声标准差。其对角为σ²当给定预测输入Xs时full返回n x m的全零矩阵——这正确反映了不同输入点之间白噪声不相关的性质。常被加在核上建模观测噪声cov_func ExpQuad(...) WhiteNoise(sigma)。平稳核Stationary所有平稳核继承Stationary基类pymc/gp/cov.py其共性是核值只依赖输入点之间的距离因此diag恒为 1.0。Stationary统一处理长度尺度参数构造函数签名如下Stationary(input_dim, lsNone, ls_invNone, active_dimsNone)ls长度尺度。input_dim 1时可以是标量列表/数组每个维度一个也支持 PyMC 随机变量input_dim 1时为标量。ls_inv逆长度尺度即1 / ls。ls与ls_inv必须且只能提供一个否则抛出ValueErrorpymc/gp/cov.py。在源码中ls_inv会被转换为ls 1.0 / ls_inv。距离计算方面square_dist先按1/ls缩放输入再计算平方距离并用pt.clip(sqd, 0.0, np.inf)消除数值误差导致的负值euclidean_dist在平方距离上开方并加上1e-12数值稳定项pymc/gp/cov.py。ExpQuad指数二次核又称平方指数 / RBF 核pymc/gp/cov.pyk(x, x) exp( -(x - x)² / (2ℓ²) )最常用的光滑核处处无穷可微。full_from_distance实现为pt.exp(-0.5 * r2)。它还实现了功率谱密度见后文谱密度一节。RatQuad有理二次核pymc/gp/cov.pyk(x, x) ( 1 (x - x)² / (2αℓ²) )^(-α)额外参数alphaα控制形状它可以看作无穷多个不同长度尺度的平方指数核的尺度混合α 越大越接近ExpQuad。源码注释给出了其作为平方指数核关于精度参数 λ ~ Gamma(α, αℓ²) 的混合表示。Exponential指数核pymc/gp/cov.pyk(x, x) exp( -||x - x|| / (2ℓ) )与 ExpQuad 不同它基于欧氏距离而非平方距离对应的是奥恩斯坦-乌伦贝克OU过程样本路径连续但不可导。Matern52 / Matern32 / Matern12Matérn 核族Matern52ν 5/2pymc/gp/cov.pyk(x, x) (1 √5·r/ℓ 5r²/(3ℓ²)) · exp(-√5·r/ℓ)Matern32ν 3/2pymc/gp/cov.pyk(x, x) (1 √3·r/ℓ) · exp(-√3·r/ℓ)Matern12ν 1/2pymc/gp/cov.pyk(x, x) exp(-r/ℓ)其中r为经长度尺度缩放后的欧氏距离。Matérn 族通过 ν 参数刻画函数可微性ν 越大函数越光滑Matern52二阶可微、Matern32一阶可微、Matern12即指数核的变体不可微。实践中Matern52是平衡光滑性与数值稳定性的常用选择。Cosine余弦核pymc/gp/cov.pyk(x, x) cos( 2π·||x - x|| / ℓ² )强周期振荡核适用于周期性但非平滑衰减的信号。Periodic周期核pymc/gp/cov.pyk(x, x) exp( -sin²(π|x - x|/T) / (2ℓ²) )需要额外参数period周期T。源码特别给出了一个易踩坑的说明PyMC 的系数约定与常见定义不同指数上是0.5而非常见的2因此当你想复现标准定义时初始化时需将长度尺度除以 2pymc/gp/cov.py。此外Periodic还实现了power_spectral_density_approx用于 HSGP 低秩近似的系数计算基于第一类修正贝塞尔函数I_j。非平稳核Linear线性核pymc/gp/cov.pyk(x, x) (x - c)(x - c)c为偏移中心点。该核产生的函数是输入空间的线性函数常用于趋势建模。Polynomial多项式核pymc/gp/cov.pyk(x, x) [ (x - c)(x - c) offset ]^d在Linear基础上增加次数d与常数项offset实现方式为对线性核的结果取(linear offset)^d。输入变换类核WarpedInputpymc/gp/cov.pyk(x, x) k_base(w(x), w(x))用任意 PyTensor 函数warp_func对输入做非线性扭曲再交给底层核cov_func计算。参数args用于向warp_func传递额外的标量或 PyMC 变量内部由handle_args包装见 pymc/gp/cov.py。典型应用将一维输入映射到高维特征空间从而获得非平稳行为。Gibbspymc/gp/cov.pyk(x, x) sqrt( 2ℓ(x)ℓ(x) / (ℓ²(x) ℓ²(x)) ) · exp( -(x - x)² / (ℓ²(x) ℓ²(x)) )使用随输入变化的长度尺度函数lengthscale_func构造非平稳核。源码明确标注仅在 1 维下测试过若active_dims长度大于 1 或input_dim ! 1会直接抛出NotImplementedErrorpymc/gp/cov.py。args同样用于传递额外参数。ScaledCovpymc/gp/cov.pyk(x, x) φ(x) · k_base(x, x) · φ(x)用非负的缩放函数φ(x)scaling_funcPyTensor 可调用对象对基础核cov_func做点级缩放实现对函数幅值的非平稳调制。源码中diag实现为cov_diag * φ(X)²full实现为outer(φ(X), φ(Xs)) * k_base(X, Xs)。多输出核CoregionCoregion协区域化核pymc/gp/cov.py用于内禀/线性协区域化模型ICM/LCM其协方差矩阵为B W·Wᵀ diag(κ)W形状(num_outputs, rank)的低秩矩阵决定各输出之间的相关性kappa形状(num_outputs,)的向量使各输出可独立变化B形状(num_outputs, num_outputs)的完整矩阵。约束条件(W, kappa)与B二者必须恰好提供一个并且该核要求恰好激活一个维度active_dims长度必须为 1否则抛出ValueErrorpymc/gp/cov.py。调用时输入应为整数索引输出编号源码将其cast为int32后通过B[index, index2]查表取值。结构核KronKronKronecker 积核pymc/gp/cov.py与普通乘法所有核共享同一份输入不同Kron 核先把输入按各因子的input_dim切分到各自子空间再对各子空间分别求核矩阵最后做 Kronecker 乘积Kron([cov1, cov2, ...])其整体input_dim等于各因子input_dim之和因子必须是协方差函数或其组合数组不支持。它通常配合gp.MarginalKron与gp.LatentKron使用用于输入为笛卡尔积结构的数据如网格数据可显著加速矩阵运算。测试中也使用pymc.math.kronecker来验证其正确性tests/gp/test_cov.py。核的代数组合可组合性是 PyMC GP 的灵魂协方差函数在 PyMC 中严格遵循核函数的代数规则——这一点被设计为语言级特性BaseCovariance 中的运算符重载用户可以像拼接乐高一样构造复杂核两个核相加仍是核Addpymc/gp/cov.pycov_func pm.gp.cov.ExpQuad(input_dim1, ls1.0) pm.gp.cov.Periodic(input_dim1, period0.5)两个核相乘仍是核Prodpymc/gp/cov.pycov_func pm.gp.cov.ExpQuad(input_dim1, ls1.0) * pm.gp.cov.Periodic(input_dim1, period0.5)核与标量的乘/加仍是核标量会自动被包装为Constantcov_func eta**2 * pm.gp.cov.Matern32(input_dim1, ls1.0) # 幅值缩放 cov_func pm.gp.cov.ExpQuad(input_dim1, ls1.0) 1.0 # 加常数偏置幂运算Exponentiatedpymc/gp/cov.pycov_func ** p要求底数必须是继承自Covariance的核Constant/WhiteNoise会报TypeError且指数必须是标量。底层实现上Combination基类pymc/gp/cov.py做了两件事自动推导元信息组合核的input_dim取各因子之交集所有因子必须一致否则报错active_dims取各因子激活维度的并集排序后。扁平化因子列表Add套Add、Prod套Prod会被自动展平避免嵌套过深。_merge_factors_cov还负责统一处理diagTrue时对数组/张量因子取对角np.diag/pt.diag的细节保证对角计算的正确性——这一点被测试用例反复验证tests/gp/test_cov.py。从测试可以看出这些运算的实际数值行为tests/gp/test_cov.pyExpQuad(1, 0.1) 1的非对角元为1.53940正是ExpQuad核值0.53940加上常数1ExpQuad(1, 0.1) WhiteNoise(sigma1)的对角元为2核对角 1 加噪声方差 1、非对角元保持0.53940。功率谱密度从频域理解平稳核对于平稳核PyMC 还实现了power_spectral_density(omega)Stationary为 HSGPHilbert Space Gaussian Process等谱方法提供理论基础ExpQuad的谱密度为高斯形式pymc/gp/cov.pyRatQuad的谱密度涉及第二类修正贝塞尔函数K_ν(z)并处理了z 0处的奇异性pymc/gp/cov.pyMatern32/Matern52的谱密度为有理函数形式pymc/gp/cov.py。组合核的谱密度遵循简单规则Add的谱密度等于各因子谱密度之和Prod仅当因子中协方差函数不多于一个时才能计算两个协方差函数相乘的谱没有解析形式会抛出NotImplementedError见 pymc/gp/cov.py。此外进行谱密度计算时要求求和项中所有核的active_dims完全一致pymc/gp/cov.py非平稳核则会收到明确的ValueError提示。在 GP 模型中使用协方差函数协方差函数本身只是零件真正发挥威力的是与 GP 实现的配合。PyMC 提供gp.Latent、gp.Marginal等多种实现统一的使用模式为先用均值函数与协方差函数实例化 GP 对象再调用prior、marginal_likelihood或conditional方法构造 PyMC 随机变量。Latent GP 的完整流程docs/source/guides/Gaussian_Processes.rstimport pymc as pm # 1. 定义核与均值 cov_func pm.gp.cov.ExpQuad(input_dim1, ls1.0) mean_func pm.gp.mean.Zero() # 2. 实例化 GP gp pm.gp.Latent(mean_func, cov_func) with pm.Model() as model: # 3. 在观测点 X 上建立函数先验 f f gp.prior(f, X) # 4. 预测点的条件分布 f* f_star gp.conditional(f_star, X_star)prior的第一个参数是随机变量名第二个是函数输入X通常来自数据也可以是 PyMC 随机变量——若输入是张量/随机变量则必须显式给出shape。gp.Marginal类没有prior方法改用marginal_likelihood需要额外提供观测数据与噪声等参数。可加 GPAdditive GP是 PyMC 的另一大特性docs/source/guides/Gaussian_Processes.rstGP 对象之间可以直接相加从而把复杂函数分解为多个独立成分并分别推断with pm.Model() as model: gp1 pm.gp.Marginal(mean_func1, cov_func1) # 长程趋势 gp2 pm.gp.Marginal(mean_func2, cov_func2) # 周期成分 gp gp1 gp2 # f1 f2 f gp.marginal_likelihood(f, X, y, noise) idata pm.sample(1000)注意相加的两个 GP 对象类型必须一致gp.Marginal不能与gp.Latent相加。若想得到某个成分如gp2的条件分布需要在conditional中以given{X: X, y: y, noise: noise, gp: gp}字典形式显式传入缓存的参数docs/source/guides/Gaussian_Processes.rstwith model: f2_star gp2.conditional(f2_star, X_star, given{X: X, y: y, noise: noise, gp: gp}) f_star gp.conditional(f_star, X_star) # 整体条件分布无需 given整体gp的条件分布无需given因为参数已在marginal_likelihood调用时被缓存而gp1、gp2自身从未调用过marginal_likelihood必须显式提供。测试驱动的工程细节与注意事项仓库中的 tests/gp/test_cov.py共 948 行是对pymc.gp.cov行为的最权威验证值得注意的工程细节包括对角一致性所有测试都校验cov(X, diagTrue)与np.diag(cov(X))严格一致如 tests/gp/test_cov.py这是BaseCovariance接口设计正确性的基本保证。标量/矩阵的左右结合a cov、cov a、M cov、cov M均被支持通过__radd__/__rmul__与__array_wrap__但形状非法的组合如三维数组会抛出ValueError: cannot combine...tests/gp/test_cov.py。组合核的元信息合并不同input_dim的核相加会报错active_dims取并集这保证组合核始终能正确切片输入。实战中的几条核心建议长度尺度选择ls可设为标量各向同性或按维度设置的数组自动相关性确定ARD当维度较多时优先用ls_inv配合正约束先验如HalfNormal、Gamma进行推断。噪声与核分离观测噪声建议显式使用WhiteNoise或gp.Marginal的noise参数保持核本身干净便于解释。周期数据优先选用Periodic注意其 0.5 缩放约定若底层核需要更强的可微性可用WrappedPeriodic将任意Stationary核卷成周期核pymc/gp/cov.py。大规模数据网格化/可分解结构数据用Kron配合gp.MarginalKron/gp.LatentKron否则可结合Periodic.power_spectral_density_approx与 HSGP 近似见 pymc/gp/hsgp_approx.py降低计算负担。小结pymc.gp.cov为贝叶斯高斯过程建模提供了完整、可组合、经过充分测试的协方差函数工具箱从经典的ExpQuad、Matérn 族、周期核到非平稳的WarpedInput、Gibbs、ScaledCov再到多输出Coregion与结构化的Kron全部遵循统一的BaseCovariance调用契约并支持加、乘、标量缩放与幂运算的代数组合。配合input_dim/active_dims的维度控制与ls/ls_inv的长度尺度参数化你可以在 pymc/gp/gp.py 提供的Latent、Marginal、MarginalKron等实现中自由搭建从简单回归到可加成分模型的各种 GP 模型。【免费下载链接】pymcBayesian Modeling and Probabilistic Programming in Python项目地址: https://gitcode.com/GitHub_Trending/py/pymc创作声明：本文部分内容由AI辅助生成（AIGC），仅供参考