非线性最小二乘与高斯-牛顿算法
开始做这个东西还是因为学校里的一次课程设计任务,找遍了全网好像也没有特别好用的,于是就自己写了一个。仅供参考。
首先,介绍下非线性最小二乘问题。
非线性最小二乘是数值最优化领域最常见的一个子问题。与一般的优化问题不同的是,非线性最小二乘的目标函数具有明确的物理意义—残差。
应用背景:消除误差或忽略无关细节,从干扰数据中提取信号或找出趋势,将大量数据降低到可管理的数量或用简单的近似来代替复杂函数。事实上,在许多时候它也不用很精确。但尽管如此,我们还是希望它能保持对原始数据的相似之处,完成这个工作的最普遍和最便于计算的方法之一就是最小二乘法。
最小二乘法是从误差拟合角度对回归模型进行参数估计或系统辨识,并在参数估计、系统辨识以及预测、预报等众多领域中得到极为广泛的应用。当我们用一个模型来描述现实中的一系列数据时,模型的预测结果与实际的测量结果总会存在一定偏差,这一偏差就称为残差。非线性最小二乘的目的就是,调整模型的参数,使得总的残差最小。
问题的分析:
若为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);