高斯过程回归（GPR）入门教程：从概念到实践

高斯过程回归（GPR）入门教程：从概念到实践

一、背景溯源：为什么需要高斯过程回归？

在机器学习中，回归问题的核心是“从数据中学习输入xxx到输出yyy的映射关系f(x)f(x)f(x)”。传统方法（如线性回归、多项式回归）的局限性很明显：

1. 需要预先假设函数形式：比如线性回归假设f(x)=wx+bf(x) = wx + bf(x)=wx+b，多项式回归假设f(x)=a0+a1x+a2x2+…f(x) = a_0 + a_1x + a_2x^2 + \dotsf(x)=a0​+a1​x+a2​x2+…，但现实中很多数据的规律是非线性、未知的，强行假设形式会导致欠拟合或过拟合；
2. 无法量化不确定性：传统模型只能给出单点预测值（如y∗=f(x∗)y^* = f(x^*)y∗=f(x∗)），但无法告诉我们“这个预测有多可靠”——而在医疗、金融、机器人等领域，不确定性估计往往比预测值更重要。

高斯过程回归（Gaussian Process Regression, GPR）作为非参数贝叶斯模型，完美解决了这两个问题：

• 它不需要假设函数形式，而是直接从数据中学习函数的分布；
• 它能输出预测值的概率分布（均值+标准差），自然量化不确定性；
• 它是贝叶斯框架下的模型，能灵活结合先验知识（如“函数是光滑的”“数据有周期性”）。

GPR的理论基础是高斯过程（Gaussian Process, GP），这一概念由统计学家Geoffrey Box和George Jenkins在1970年代提出，后来被机器学习社区（如Carl Edward Rasmussen）推广到回归任务中，成为小样本、高价值数据场景的“瑞士军刀”。

二、核心思想：把“函数”当作“随机变量”

要理解GPR，首先要打破“函数是固定公式”的思维定式——GPR把所有可能的函数 $ f(x) $ 看作一个“函数空间”里的随机变量。具体来说：

• 对于任意一个输入点xxx，函数值f(x)f(x)f(x)是一个随机变量；
• 对于任意kkk个输入点x1,x2,…,xkx_1, x_2, \dots, x_kx1​,x2​,…,xk​，对应的函数值f(x1),f(x2),…,f(xk)f(x_1), f(x_2), \dots, f(x_k)f(x1​),f(x2​),…,f(xk​)的联合分布是k维高斯分布。

这种“任意有限点的联合分布都是高斯”的性质，就是高斯过程的定义。换句话说：

高斯过程是“函数的分布”，由均值函数和协方差函数完全描述。

1. 高斯过程的数学定义

一个高斯过程记为：f(x)∼GP(m(x),k(x,x′)) f(x) \sim \mathcal{GP}\left( m(x), k(x, x') \right)f(x)∼GP(m(x),k(x,x′))其中：

• 均值函数m(x)=E[f(x)]m(x) = \mathbb{E}\left[ f(x) \right]m(x)=E[f(x)]：描述所有可能函数在点 $ x $ 处的平均取值（比如“大多数函数在 $ x=0 $ 处的值接近0”）；
• 协方差函数k(x,x′)=Cov[f(x),f(x′)]=E[(f(x)−m(x))(f(x′)−m(x′))]k(x, x') = \text{Cov}\left[ f(x), f(x') \right] = \mathbb{E}\left[ (f(x)-m(x))(f(x')-m(x')) \right]k(x,x′)=Cov[f(x),f(x′)]=E[(f(x)−m(x))(f(x′)−m(x′))]：描述两个点xxx和x′x'x′处函数值的相关性（比如“相邻点的函数值更接近”）。

2. 关键类比：从“有限维高斯”到“无限维高斯”

我们已经熟悉有限维高斯分布（比如一维的N(μ,σ2)N(\mu, \sigma^2)N(μ,σ2)，二维的N(μ,Σ)N(\boldsymbol{\mu}, \boldsymbol{\Sigma})N(μ,Σ)），而高斯过程是无限维的高斯分布——它覆盖了所有可能的输入点xxx，对应的函数值f(x)f(x)f(x)构成一个无限长的随机向量。

但高斯过程的“无限性”并不难处理，因为我们只需要关注有限个点的联合分布：比如训练集的NNN个点X={x1,…,xN}X = \{x_1, \dots, x_N\}X={x1​,…,xN​}，对应的函数值f(X)=[f(x1),…,f(xN)]Tf(X) = [f(x_1), \dots, f(x_N)]^Tf(X)=[f(x1​),…,f(xN​)]T的分布是：f(X)∼N(m(X),K(X,X)) f(X) \sim \mathcal{N}\left( m(X), K(X, X) \right)f(X)∼N(m(X),K(X,X))其中：

• m(X)=[m(x1),…,m(xN)]Tm(X) = [m(x_1), \dots, m(x_N)]^Tm(X)=[m(x1​),…,m(xN​)]T是均值向量；
• K(X,X)K(X, X)K(X,X)是N×NN \times NN×N的协方差矩阵（也叫Gram矩阵），元素为Ki,j=k(xi,xj)K_{i,j} = k(x_i, x_j)Ki,j​=k(xi​,xj​)。

三、算法原理：从先验到后验的贝叶斯更新

GPR的核心流程是贝叶斯推断：用训练数据更新对函数f(x)f(x)f(x)的先验分布，得到后验分布，再用后验分布预测新点。

1. 步骤1：定义先验分布

先验分布是我们对函数f(x)f(x)f(x)的初始假设，由均值函数和协方差函数决定。常见选择：

• 均值函数：通常设为0（m(x)=0m(x) = 0m(x)=0），因为数据可以通过中心化（减去均值）消除偏移；
• 协方差函数：选择能刻画数据规律的形式（如光滑性、周期性），最常用的是平方指数（Squared Exponential, SE）协方差函数：kSE(x,x′)=σf2exp⁡(−∥x−x′∥22l2) k_{\text{SE}}(x, x') = \sigma_f^2 \exp\left( -\frac{\|x - x'\|^2}{2l^2} \right)kSE​(x,x′)=σf2​exp(−2l2∥x−x′∥2​)其中：
  • σf2\sigma_f^2σf2​：信号方差，控制函数的波动幅度（越大，函数越“起伏”）；
  • lll：长度尺度，控制函数随 $ x $ 变化的快慢（越小，函数变化越“剧烈”）；
  • ∥x−x′∥2\|x - x'\|^2∥x−x′∥2：输入点 $ x $ 和 $ x’ $ 的欧氏距离平方。

SE协方差函数的特点是光滑性：相邻点的函数值高度相关，因此先验中的函数都是光滑曲线。

2. 步骤2：加入观测噪声

训练数据中的输出yiy_iyi​不是纯函数值f(xi)f(x_i)f(xi​)，而是带有噪声：yi=f(xi)+ϵi y_i = f(x_i) + \epsilon_iyi​=f(xi​)+ϵi​其中ϵi∼N(0,σn2)\epsilon_i \sim \mathcal{N}(0, \sigma_n^2)ϵi​∼N(0,σn2​)是独立同分布的高斯噪声（σn2\sigma_n^2σn2​是噪声方差）。

因此，训练集的输出向量y=[y1,…,yN]T\boldsymbol{y} = [y_1, \dots, y_N]^Ty=[y1​,…,yN​]T的分布是：y∼N(m(X),K(X,X)+σn2I) \boldsymbol{y} \sim \mathcal{N}\left( m(X), K(X, X) + \sigma_n^2 I \right)y∼N(m(X),K(X,X)+σn2​I)其中III是N×NN \times NN×N单位矩阵（噪声独立，协方差矩阵对角线为σn2\sigma_n^2σn2​）。

3. 步骤3：计算后验分布（核心推导）

根据贝叶斯定理，后验分布p(f∣X,y)p(f | X, \boldsymbol{y})p(f∣X,y)满足：p(f∣X,y)∝p(y∣f,X)⋅p(f∣X) p(f | X, \boldsymbol{y}) \propto p(\boldsymbol{y} | f, X) \cdot p(f | X)p(f∣X,y)∝p(y∣f,X)⋅p(f∣X)其中：

• p(f∣X)p(f | X)p(f∣X)：先验分布（高斯）；
• p(y∣f,X)p(\boldsymbol{y} | f, X)p(y∣f,X)：似然（高斯，因为噪声是高斯的）。

关键性质：高斯分布的共轭性——先验和似然都是高斯，后验也一定是高斯！

我们需要预测新点x∗x^*x∗的函数值f(x∗)f(x^*)f(x∗)，因此考虑联合分布：[yf(x∗)]∼N([m(X)m(x∗)],[K(X,X)+σn2IK(X,x∗)K(x∗,X)k(x∗,x∗)]) \begin{bmatrix} \boldsymbol{y} \\ f(x^*) \end{bmatrix} \sim \mathcal{N}\left( \begin{bmatrix} m(X) \\ m(x^*) \end{bmatrix}, \begin{bmatrix} K(X,X) + \sigma_n^2 I & K(X, x^*) \\ K(x^*, X) & k(x^*, x^*) \end{bmatrix} \right)[yf(x∗)​]∼N([m(X)m(x∗)​],[K(X,X)+σn2​IK(x∗,X)​K(X,x∗)k(x∗,x∗)​])其中：

• K(X,x∗)=[k(x1,x∗),…,k(xN,x∗)]TK(X, x^*) = [k(x_1, x^*), \dots, k(x_N, x^*)]^TK(X,x∗)=[k(x1​,x∗),…,k(xN​,x∗)]T：训练点与新点的协方差向量；
• K(x∗,X)=K(X,x∗)TK(x^*, X) = K(X, x^*)^TK(x∗,X)=K(X,x∗)T：转置；
• k(x∗,x∗)k(x^*, x^*)k(x∗,x∗)：新点自身的协方差（即σf2\sigma_f^2σf2​，因为x∗=x′x^* = x'x∗=x′）。

根据高斯分布的条件概率公式：若联合分布[Z1Z2]∼N([μ1μ2],[Σ11Σ12Σ21Σ22])\begin{bmatrix} Z_1 \\ Z_2 \end{bmatrix} \sim \mathcal{N}\left( \begin{bmatrix} \mu_1 \\ \mu_2 \end{bmatrix}, \begin{bmatrix} \Sigma_{11} & \Sigma_{12} \\ \Sigma_{21} & \Sigma_{22} \end{bmatrix} \right)[Z1​Z2​​]∼N([μ1​μ2​​],[Σ11​Σ21​​Σ12​Σ22​​])，则条件分布Z2∣Z1=z1Z_2 | Z_1 = z_1Z2​∣Z1​=z1​为：Z2∣Z1=z1∼N(μ2+Σ21Σ11−1(z1−μ1),Σ22−Σ21Σ11−1Σ12) Z_2 | Z_1 = z_1 \sim \mathcal{N}\left( \mu_2 + \Sigma_{21} \Sigma_{11}^{-1} (z_1 - \mu_1), \Sigma_{22} - \Sigma_{21} \Sigma_{11}^{-1} \Sigma_{12} \right)Z2​∣Z1​=z1​∼N(μ2​+Σ21​Σ11−1​(z1​−μ1​),Σ22​−Σ21​Σ11−1​Σ12​)

将此公式应用到我们的联合分布中（令Z1=y,Z2=f(x∗)Z_1 = \boldsymbol{y}, Z_2 = f(x^*)Z1​=y,Z2​=f(x∗)），可得到后验分布的均值和协方差：

（1）后验均值（预测值）

mpost(x∗)=m(x∗)+K(x∗,X)(K(X,X)+σn2I)−1(y−m(X)) m_{\text{post}}(x^*) = m(x^*) + K(x^*, X) \left( K(X,X) + \sigma_n^2 I \right)^{-1} \left( \boldsymbol{y} - m(X) \right)mpost​(x∗)=m(x∗)+K(x∗,X)(K(X,X)+σn2​I)−1(y−m(X))

（2）后验协方差（不确定性）

kpost(x∗,x∗)=k(x∗,x∗)−K(x∗,X)(K(X,X)+σn2I)−1K(X,x∗) k_{\text{post}}(x^*, x^*) = k(x^*, x^*) - K(x^*, X) \left( K(X,X) + \sigma_n^2 I \right)^{-1} K(X, x^*)kpost​(x∗,x∗)=k(x∗,x∗)−K(x∗,X)(K(X,X)+σn2​I)−1K(X,x∗)

4. 步骤4：预测输出（含噪声）

如果要预测新点的输出y∗y^*y∗（而非纯函数值f(x∗)f(x^*)f(x∗)），需要加上新的噪声ϵ∗∼N(0,σn2)\epsilon^* \sim \mathcal{N}(0, \sigma_n^2)ϵ∗∼N(0,σn2​)，因此：y∗∣X,y,x∗∼N(mpost(x∗),kpost(x∗,x∗)+σn2) y^* | X, \boldsymbol{y}, x^* \sim \mathcal{N}\left( m_{\text{post}}(x^*), k_{\text{post}}(x^*, x^*) + \sigma_n^2 \right)y∗∣X,y,x∗∼N(mpost​(x∗),kpost​(x∗,x∗)+σn2​)

解释：

• 预测值是后验均值mpost(x∗)m_{\text{post}}(x^*)mpost​(x∗)；
• 预测的不确定性是后验标准差kpost(x∗,x∗)+σn2\sqrt{k_{\text{post}}(x^*, x^*) + \sigma_n^2}kpost​(x∗,x∗)+σn2​​——这个值越大，预测越不可靠。

四、完整模型求解步骤（小白可落地）

GPR的求解过程可拆分为7个具体步骤，每一步都有明确的操作指南：

步骤1：数据准备

• 收集训练数据D={(x1,y1),(x2,y2),…,(xN,yN)}D = \{(x_1, y_1), (x_2, y_2), \dots, (x_N, y_N)\}D={(x1​,y1​),(x2​,y2​),…,(xN​,yN​)}，其中xix_ixi​是输入特征（可以是一维或多维），yiy_iyi​是输出标签；
• 对输入特征进行标准化（如xi′=(xi−μx)/σxx_i' = (x_i - \mu_x)/\sigma_xxi′​=(xi​−μx​)/σx​）：因为协方差函数（如SE）对特征尺度敏感，标准化后能加速超参数优化；
• 对输出标签进行中心化（如yi′=yi−μyy_i' = y_i - \mu_yyi′​=yi​−μy​）：方便后续设均值函数m(x)=0m(x) = 0m(x)=0（预测时再还原均值）。

步骤2：选择均值函数与协方差函数

• 均值函数：优先选m(x)=0m(x) = 0m(x)=0（若数据已中心化），或简单线性函数（如m(x)=a+bxm(x) = a + bxm(x)=a+bx）；
• 协方差函数：根据数据规律选择（以下为常见类型）：

步骤3：初始化超参数

协方差函数和噪声的参数称为超参数（如SE的 $ \sigma_f^2, l $，噪声的 $ \sigma_n^2 $）。初始值可根据数据经验设置：

• σf2\sigma_f^2σf2​：设为输出数据的方差（$ \sigma_f^2 = \text{Var}(y) $）；
• lll：设为输入特征的标准差（$ l = \sigma_x $）；
• σn2\sigma_n^2σn2​：设为小正数（如 $ 0.1 $），避免过拟合。

步骤4：优化超参数（关键！）

超参数直接决定模型性能，需通过最大化边际似然（Marginal Likelihood）优化。边际似然是对函数fff积分后的似然（消除了fff的影响）：p(y∣X,θ)=∫p(y∣f,X,θ)p(f∣X,θ)df p(\boldsymbol{y} | X, \theta) = \int p(\boldsymbol{y} | f, X, \theta) p(f | X, \theta) dfp(y∣X,θ)=∫p(y∣f,X,θ)p(f∣X,θ)df其中θ\thetaθ是超参数集合（如 $ \theta = {\sigma_f^2, l, \sigma_n^2} $）。

对高斯分布，边际似然的对数形式（更易计算）为：log⁡p(y∣X,θ)=−12(y−m(X))TΣ−1(y−m(X))−12log⁡∣Σ∣−N2log⁡(2π) \log p(\boldsymbol{y} | X, \theta) = -\frac{1}{2} (\boldsymbol{y} - m(X))^T \Sigma^{-1} (\boldsymbol{y} - m(X)) - \frac{1}{2} \log|\Sigma| - \frac{N}{2} \log(2\pi)logp(y∣X,θ)=−21​(y−m(X))TΣ−1(y−m(X))−21​log∣Σ∣−2N​log(2π)其中Σ=K(X,X;θ)+σn2I\Sigma = K(X,X; \theta) + \sigma_n^2 IΣ=K(X,X;θ)+σn2​I（训练集的协方差矩阵+噪声）。

优化方法：

• 用数值优化算法（如L-BFGS、梯度下降）最大化log⁡p(y∣X,θ)\log p(\boldsymbol{y} | X, \theta)logp(y∣X,θ)；
• 计算对数边际似然对每个超参数的梯度（需用矩阵求导规则），指导优化方向。

步骤5：计算后验均值与协方差

用优化后的超参数θ∗\theta^*θ∗，计算：
3. 训练集协方差矩阵：K(X,X;θ∗)K(X,X; \theta^*)K(X,X;θ∗)；
4. 新点与训练集的协方差向量：K(x∗,X;θ∗)K(x^*, X; \theta^*)K(x∗,X;θ∗)；
5. 后验均值（预测值）：mpost(x∗)=K(x∗,X;θ∗)Σ−1ym_{\text{post}}(x^*) = K(x^*, X; \theta^*) \Sigma^{-1} \boldsymbol{y}mpost​(x∗)=K(x∗,X;θ∗)Σ−1y（若m(x)=0m(x)=0m(x)=0）；
6. 后验协方差（不确定性）：kpost(x∗,x∗)=k(x∗,x∗;θ∗)−K(x∗,X;θ∗)Σ−1K(X,x∗;θ∗)k_{\text{post}}(x^*, x^*) = k(x^*, x^*; \theta^*) - K(x^*, X; \theta^*) \Sigma^{-1} K(X, x^*; \theta^*)kpost​(x∗,x∗)=k(x∗,x∗;θ∗)−K(x∗,X;θ∗)Σ−1K(X,x∗;θ∗)。

步骤6：预测新点

对每个新输入x∗x^*x∗，输出：

• 点预测：y^∗=mpost(x∗)+μy\hat{y}^* = m_{\text{post}}(x^*) + \mu_yy^​∗=mpost​(x∗)+μy​（还原中心化的均值）；
• 不确定性区间：y^∗±2kpost(x∗,x∗)+σn2\hat{y}^* \pm 2\sqrt{k_{\text{post}}(x^*, x^*) + \sigma_n^2}y^​∗±2kpost​(x∗,x∗)+σn2​​（95%置信区间）。

步骤7：评估模型性能

用验证集（如将训练集按8:2拆分）评估模型：

• 回归指标：均方误差（MSE）、平均绝对误差（MAE）；
• 不确定性指标：负对数概率（NLL）——衡量预测分布与真实值的匹配程度，越小越好：NLL=−1M∑i=1Mlog⁡p(yi∗∣X,y,xi∗) \text{NLL} = -\frac{1}{M} \sum_{i=1}^M \log p(y_i^* | X, \boldsymbol{y}, x_i^*)NLL=−M1​i=1∑M​logp(yi∗​∣X,y,xi∗​)其中 $ M $ 是验证集大小。

五、适用边界：GPR的“能”与“不能”

GPR是强大的工具，但并非万能，需明确其适用场景和局限性。

1. 优势（适用场景）

• 小数据集：计算复杂度为O(N3)O(N^3)O(N3)（协方差矩阵求逆），小样本（N<1000N < 1000N<1000）时快速且准确；
• 需要不确定性估计：如医疗诊断（“预测患者康复率为70%，置信区间60%-80%”）、金融风险（“股票涨幅预测的标准差为2%”）；
• 非线性关系：无需假设函数形式，能拟合复杂曲线（如正弦函数、指数函数）；
• 结合先验知识：通过协方差函数融入领域知识（如“数据有周期性”“相邻点相关”）。

2. 局限性（不适用场景）

• 大数据集：当N>104N > 10^4N>104时，O(N3)O(N^3)O(N3)的计算量会导致内存溢出或计算超时；
• 高维数据：输入特征维度d>20d > 20d>20时，“维度灾难”会让协方差函数无法捕捉相关性（高维空间中所有点的距离都很大）；
• 非高斯噪声：GPR假设噪声是高斯的，若噪声是重尾（如 outliers 多）或非对称的，性能会下降；
• 超参数敏感：协方差函数选择或超参数优化不当，会导致过拟合（如lll太小，函数波动过大）或欠拟合（如lll太大，函数太光滑）。

六、总结：GPR的核心逻辑

高斯过程回归的本质是用高斯分布描述函数的不确定性，通过贝叶斯更新从先验到后验，最终输出“预测值+不确定性”。

对小白来说，最关键的是理解：

• 高斯过程是“函数的分布”，由均值和协方差函数定义；
• 后验分布是条件高斯分布，均值是预测值，协方差是不确定性；
• 协方差函数和超参数决定了模型的“先验假设”，是性能的关键。

通过一个简单的一维曲线拟合例子（如用SE协方差函数拟合y=sin⁡(x)+噪声y = \sin(x) + \text{噪声}y=sin(x)+噪声），你可以直观看到GPR的效果：先验是一组光滑曲线，训练数据后，后验曲线会“收缩”到训练点附近，且预测点的不确定性随远离训练数据而增大——这正是GPR的魅力所在！

案例解析

importnumpyasnpimportmatplotlib.pyplotaspltfromscipy.optimizeimportminimize# ------------------------------# 数据生成函数：生成带噪声的正弦曲线数据# ------------------------------defgenerate_data(n_samples=50,noise_std=0.1):# 生成输入x：在0到4π之间均匀采样X=np.linspace(0,4*np.pi,n_samples)[:,np.newaxis]# 生成真实函数值y = sin(x)y_true=np.sin(X).flatten()# 生成带高斯噪声的观测值y_noisy=y_true+np.random.normal(0,noise_std,size=n_samples)returnX,y_noisy,y_true# ------------------------------# 平方指数(SE)协方差函数# ------------------------------defse_covariance(X1,X2,sigma_f,length_scale):# 计算欧氏距离的平方distance=np.sum(X1**2,1)[:,np.newaxis]+np.sum(X2**2,1)-2*np.dot(X1,X2.T)# 平方指数协方差矩阵K=sigma_f**2*np.exp(-0.5*distance/(length_scale**2))returnK# ------------------------------# 对数边际似然函数（用于超参数优化）# ------------------------------deflog_marginal_likelihood(theta,X,y,sigma_n):# 解析超参数：sigma_f和length_scalesigma_f=np.exp(theta[0])length_scale=np.exp(theta[1])# 计算训练集协方差矩阵K + σn²IK=se_covariance(X,X,sigma_f,length_scale)+sigma_n**2*np.eye(X.shape[0])# 计算K的逆和行列式K_inv=np.linalg.inv(K)log_det_K=np.linalg.slogdet(K)[1]# 对数边际似然公式lml=-0.5*y.T @ K_inv @ y-0.5*log_det_K-0.5*X.shape[0]*np.log(2*np.pi)# 返回负的lml，因为优化器默认最小化return-lml# ------------------------------# 高斯过程回归预测函数# ------------------------------defgpr_predict(X_train,y_train,X_test,theta_opt,sigma_n):# 解析优化后的超参数sigma_f_opt=np.exp(theta_opt[0])length_scale_opt=np.exp(theta_opt[1])# 计算训练集协方差矩阵K_train=se_covariance(X_train,X_train,sigma_f_opt,length_scale_opt)+sigma_n**2*np.eye(X_train.shape[0])# 计算测试集与训练集的协方差向量K_test_train=se_covariance(X_test,X_train,sigma_f_opt,length_scale_opt)# 计算测试集自身的协方差矩阵K_test=se_covariance(X_test,X_test,sigma_f_opt,length_scale_opt)# 计算后验均值和协方差K_train_inv=np.linalg.inv(K_train)mu_post=K_test_train @ K_train_inv @ y_train cov_post=K_test-K_test_train @ K_train_inv @ K_test_train.T# 预测值的方差（对角线元素）var_post=np.diag(cov_post)# 带噪声的预测方差（用于置信区间）var_post_noisy=var_post+sigma_n**2returnmu_post,np.sqrt(var_post_noisy)# ------------------------------# 主程序：数据生成、模型训练与预测# ------------------------------if__name__=="__main__":# 步骤1：生成训练数据X_train,y_train,y_true_train=generate_data(n_samples=50,noise_std=0.1)# 步骤2：数据预处理# 输入标准化X_mean=np.mean(X_train)X_std=np.std(X_train)X_train_scaled=(X_train-X_mean)/X_std# 输出中心化y_mean=np.mean(y_train)y_train_centered=y_train-y_mean# 步骤3：初始化超参数（取对数后的值，便于优化）initial_theta=[np.log(np.std(y_train)),np.log(np.std(X_train))]# [log(sigma_f), log(length_scale)]sigma_n=0.1# 噪声方差（初始值）# 步骤4：优化超参数result=minimize(fun=lambdatheta:log_marginal_likelihood(theta,X_train_scaled,y_train_centered,sigma_n),x0=initial_theta,method='L-BFGS-B')theta_opt=result.x# 步骤5：生成测试数据（用于预测）X_test=np.linspace(0,4*np.pi,200)[:,np.newaxis]X_test_scaled=(X_test-X_mean)/X_std# 步骤6：进行预测mu_post_centered,std_post=gpr_predict(X_train_scaled,y_train_centered,X_test_scaled,theta_opt,sigma_n)# 还原中心化的预测值mu_post=mu_post_centered+y_mean# 步骤7：可视化结果plt.figure(figsize=(10,6))# 真实函数曲线plt.plot(X_test,np.sin(X_test).flatten(),'r--',label='True Function (sin(x))')# 训练数据点plt.scatter(X_train,y_train,color='b',label='Training Data (Noisy)')# GPR预测均值plt.plot(X_test,mu_post,'k-',label='GPR Prediction Mean')# 95%置信区间（均值±2倍标准差）plt.fill_between(X_test.flatten(),mu_post-2*std_post,mu_post+2*std_post,alpha=0.1,color='k',label='95% Confidence Interval')plt.xlabel('Input x')plt.ylabel('Output y')plt.title('Gaussian Process Regression - Sinusoidal Curve Fitting')plt.legend()plt.grid(True)plt.show()

OK，今天数学建模算法分享就到这里结束了，若内容对大家有所帮助,可以收藏慢慢看,感谢大家支持
本文有若有不足之处，希望各位帅哥美女们能给出宝贵的意见。谢谢大家！
初来乍到，本期制作不易希望各位帅哥美女们能动动小手，三连走一走！！！
如果大家觉得主播写的还不错的话可以关注一下，每日都会有推文的，期待一下明天的干货把，嘿嘿~