非线性最小二乘问题的高斯-牛顿算法

@非线性最小二乘问题的高斯-牛顿算法

非线性最小二乘与高斯-牛顿算法

开始做这个东西还是因为学校里的一次课程设计任务，找遍了全网好像也没有特别好用的，于是就自己写了一个。仅供参考。

首先，介绍下非线性最小二乘问题。

非线性最小二乘是数值最优化领域最常见的一个子问题。与一般的优化问题不同的是，非线性最小二乘的目标函数具有明确的物理意义—残差。

应用背景：消除误差或忽略无关细节，从干扰数据中提取信号或找出趋势，将大量数据降低到可管理的数量或用简单的近似来代替复杂函数。事实上，在许多时候它也不用很精确。但尽管如此，我们还是希望它能保持对原始数据的相似之处，完成这个工作的最普遍和最便于计算的方法之一就是最小二乘法。
最小二乘法是从误差拟合角度对回归模型进行参数估计或系统辨识，并在参数估计、系统辨识以及预测、预报等众多领域中得到极为广泛的应用。当我们用一个模型来描述现实中的一系列数据时，模型的预测结果与实际的测量结果总会存在一定偏差，这一偏差就称为残差。非线性最小二乘的目的就是，调整模型的参数，使得总的残差最小。

问题的分析：
若为x的非线性函数，非线性导致无法直接写出导数形式以及准确求得函数的全局最优解，因此只能通过求目标函数的局部最小值，并且尝试从局部最小值拓展到全局最小值。局部最小值的求解过程中需要涉及到驻点以及目标函数的一、二阶导数，故可以将该问题转化为多元函数的极值问题。

对于非线性最小二乘问题，一般有LM法、信赖域方法、梯度下降法、高斯牛顿法等处理方法。在本文中，主要采用高斯牛顿法研究非线性最小二乘问题。

高斯-牛顿算法

高斯牛顿法对里的每一个函数进行了一阶泰勒近似，将原本的非线性函数变成了线性函数。当其线性化完成后，原本的非线性最小二乘问题也就成为了线性最小二乘。非线性函数通过一阶泰勒展开舍弃高阶项可以将其变成线性函数，即非线性函数在x处的一阶泰勒近似使其线性化。因为线性化函数是非常容易操作的，而非线性函数就非常复杂。所以，高斯牛顿法的实质就是将非线性函数通过泰勒展开近似成线性函数，然后进行求解。

MATLAB代码

以下代码仅适用于这道题。

function[X]=CHH_GaussNewton(X,K,s)
format long;
k=0;  
x=1:1:8;
p=[8.1 11.0 14.7 19.7 26.7 35.2 44.4 55.9];
while k<K
    R=r_function(X);
    J=Jacobian(X);
    judge=inv(J'*J)*J'*R
   if norm(judge,2)<s
       fprintf('经过%d次迭代：最优参数X=[%f,%f]\n',k,X(1),X(2)')
       break
   end
   X=X-judge;
   k=k+1;
   fprintf('第%d迭代：X=[%f,%f]\n',k,X(1),X(2));
  A=X(1);
B=X(2);
y=A.*exp(B.*x);
subplot(3,4,k);
plot(x,p,'*');
hold on
plot(x,y);
xlabel('time/year');
ylabel('population/million');
title(['第',num2str(k),'次迭代拟合曲线图']);
end

function[J]=Jacobian(X)
A=X(1);
B=X(2);
J=[exp(B),A.*exp(B);exp(2.*B),2.*A.*exp(2.*B);exp(3.*B),3.*A.*exp(3.*B);exp(4.*B),4.*A.*exp(4.*B);exp(5.*B),5.*A.*exp(5.*B);exp(6.*B),6.*A.*exp(6.*B);exp(7.*B),7.*A.*exp(7.*B);exp(8.*B),8.*A.*exp(8.*B)]
``

function[r]=r_function(X)
A=X(1);
B=X(2);
r=[A.*exp(B)-8.1;A.*exp(2.*B)-11;A.*exp(3.*B)-14.7;A.*exp(4.*B)-19.7;A.*exp(5.*B)-26.7;A.*exp(6.*B)-35.2;A.*exp(7.*B)-44.4;A.*exp(8.*B)-55.9];

clc;
clear;
format long;
s=1e-3;
X=[1,1];
K=100;
CHH_GaussNewton(X,K,s);