考虑以下曲线(Rat43) 的拟合问题:
$$y=\frac{b_1}{(1+e^{b_2-b_{3}x}{)}^{1/b_4}}$$
也就是说，给定一些数据 { x i , y i } , ∀ i = 1 , . . . , n \{x_i, y_i\},\ \forall i=1,... ,n {xi​,yi​}, ∀i=1,...,n，确定最适合该数据的参数$b_1,b_2,b_3,b_4$。

我们面临的问题就是求解 $b_1,b_2,b_3,b_4$ 使下列表达式的取值最小：
$$\begin{array}{rl}E(b_1,b_2,b_3,b_4)& =\sum _i{f}^2(b_1,b_2,b_3,b_4;x_i,y_i)\\ & =\sum _i{\left(\frac{b_1}{(1+e^{b_2-b_3{x}_i}{)}^{1/b_4}}-y_i\right)}^2\end{array}$$

最佳拟合的概念取决于用来衡量拟合质量的目标函数的选择，而该目标函数又取决于产生观测结果的潜在噪声过程。当噪声是高斯噪声时，将差的平方和最小化是正确的做法。在这种情况下，参数的最优值是最大似然估计。

为了使用Ceres求解器来解决这个问题，我们需要定义一个CostFunction来计算给定x和y的残差f及其对b1,b2,b3和b4的导数。根据高等数学中微积分的知识，我们可以算出f的一系列导数：
$$\begin{array}{rl}{D}_{1}f(b_1,b_2,b_3,b_4;x,y)& =\frac{1}{(1+e^{b_2-b_{3}x}{)}^{1/b_4}}\\ D_{2}f(b_1,b_2,b_3,b_4;x,y)& =\frac{-b_1{e}^{b_2-b_{3}x}}{b_4(1+e^{b_2-b_{3}x}{)}^{1/b_4+1}}\\ D_{3}f(b_1,b_2,b_3,b_4;x,y)& =\frac{b_{1}x{e}^{b_2-b_{3}x}}{b_4(1+e^{b_2-b_{3}x}{)}^{1/b_4+1}}\\ D_{4}f(b_1,b_2,b_3,b_4;x,y)& =\frac{b_1\log \left(1+e^{b_2-b_{3}x}\right)}{b_4^2(1+e^{b_2-b_{3}x}{)}^{1/b_4}}\end{array}$$

根据这些手动计算出的导数，我们现在可以实现CostFunction了:

class Rat43Analytic : public SizedCostFunction<1,4> {public:Rat43Analytic(const double x, const double y) : x_(x), y_(y) {}virtual ~Rat43Analytic() {}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) return true;double* jacobian = jacobians[0];if (!jacobian) return true;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 ~Rat43AnalyticOptimized() {}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) return true;double* jacobian = jacobians[0];if (!jacobian) return true;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_;};

这两种实现在性能上有什么不同?

Rat43AnalyticOptimized比Rat43Analytic快2.8倍。这种运行时上的差异并不少见。为了从分析计算的导数中获得最佳性能，通常需要优化代码以考虑公共子表达式。

什么时候应该使用analytical derivatives?

• 1.表达式很简单，例如大部分是线性的

• 2.像Maple、Mathematica或symy这样的计算机代数系统可以用来符号化地区分目标函数，并生成c++来计算它们。

• 3.式子中有一些代数结构可以实现比自动微分有更好的性能。
  也就是说, 获得在计算倒数之外的最大性能需要大量的工作.在沿着这条路径走下去之前，评估雅可比矩阵的计算花费是整个求解时间的一小部分是很有用的，，记住Amdahl法则是你的朋友。

• 4.没有其他方法来计算导数，例如你想要计算一个多项式的根的导数:
  $$a_3(x,y)z^3+a_2(x,y)z^2+a_1(x,y)z+a_0(x,y)=0$$

  对于x，y，这就需要用到反函数定理

• 5.你喜欢链式法则也喜欢手动计算代数运算。