Ceres Solver中的解析导数实现与优化指南-优快云博客

本文链接：https://blog.youkuaiyun.com/gitblog_01077/article/details/148549136

Ceres Solver中的解析导数实现与优化指南

ceres-solver A large scale non-linear optimization library 项目地址: https://gitcode.com/gh_mirrors/ce/ceres-solver

引言

在非线性最小二乘问题求解中，Ceres Solver是一个功能强大的库。本文将深入探讨如何使用解析导数（Analytic Derivatives）来实现和优化代价函数，特别是在处理复杂数学模型时的最佳实践。

问题背景

我们以Rat43曲线拟合问题为例，该曲线模型为：

y = b₁ / (1 + e^(b₂ - b₃x))^(1/b₄)

给定一组数据点{(xᵢ, yᵢ)}，我们的目标是找到参数b₁, b₂, b₃和b₄的最佳拟合值，使得以下目标函数最小化：

E(b₁, b₂, b₃, b₄) = Σ [b₁ / (1 + e^(b₂ - b₃xᵢ))^(1/b₄) - yᵢ]²

解析导数实现

导数计算

通过微分计算，我们得到各参数的偏导数：

对b₁的偏导数：1 / (1 + e^(b₂ - b₃x))^(1/b₄)
对b₂的偏导数：-b₁e^(b₂ - b₃x) / [b₄(1 + e^(b₂ - b₃x))^(1/b₄ + 1)]
对b₃的偏导数：x·b₁e^(b₂ - b₃x) / [b₄(1 + e^(b₂ - b₃x))^(1/b₄ + 1)]
对b₄的偏导数：b₁·ln(1 + e^(b₂ - b₃x)) / [b₄²(1 + e^(b₂ - b₃x))^(1/b₄)]

基础实现

基于上述导数，我们可以实现基础的CostFunction：

class Rat43Analytic : public SizedCostFunction<1,4> {
   public:
     Rat43Analytic(const double x, const double y) : x_(x), y_(y) {}
     virtual bool Evaluate(double const* const* parameters,
                          double* residuals,
                          double** jacobians) const {
       const double b1 = parameters[0][0];
       const double b2 = parameters[0][1];
       const double b3 = parameters[0][2];
       const double b4 = parameters[0][3];

       residuals[0] = b1 * pow(1 + exp(b2 - b3 * x_), -1.0 / b4) - y_;

       if (jacobians && jacobians[0]) {
         double* jacobian = jacobians[0];
         jacobian[0] = pow(1 + exp(b2 - b3 * x_), -1.0 / b4);
         jacobian[1] = -b1 * exp(b2 - b3 * x_) *
                       pow(1 + exp(b2 - b3 * x_), -1.0 / b4 - 1) / b4;
         jacobian[2] = x_ * b1 * exp(b2 - b3 * x_) *
                       pow(1 + exp(b2 - b3 * x_), -1.0 / b4 - 1) / b4;
         jacobian[3] = b1 * log(1 + exp(b2 - b3 * x_)) *
                       pow(1 + exp(b2 - b3 * x_), -1.0 / b4) / (b4 * b4);
       }
       return true;
     }

    private:
     const double x_;
     const double y_;
};

优化实现

上述实现存在重复计算，我们可以通过缓存中间结果来优化：

class Rat43AnalyticOptimized : public SizedCostFunction<1,4> {
   public:
     Rat43AnalyticOptimized(const double x, const double y) : x_(x), y_(y) {}
     virtual bool Evaluate(double const* const* parameters,
                          double* residuals,
                          double** jacobians) const {
       const double b1 = parameters[0][0];
       const double b2 = parameters[0][1];
       const double b3 = parameters[0][2];
       const double b4 = parameters[0][3];

       const double t1 = exp(b2 - b3 * x_);
       const double t2 = 1 + t1;
       const double t3 = pow(t2, -1.0 / b4);
       residuals[0] = b1 * t3 - y_;

       if (jacobians && jacobians[0]) {
         double* jacobian = jacobians[0];
         const double t4 = pow(t2, -1.0 / b4 - 1);
         jacobian[0] = t3;
         jacobian[1] = -b1 * t1 * t4 / b4;
         jacobian[2] = -x_ * jacobian[1];
         jacobian[3] = b1 * log(t2) * t3 / (b4 * b4);
       }
       return true;
     }

   private:
     const double x_;
     const double y_;
};