前言

最小二乘线性拟合是常规操作，本文直接跨过。由于多元线性拟合和一元非线性拟合关系密切，故本文将其二放在一起讨论。本文重点是实现单元的非线性拟合。最后补充一句：不论是一元线性、一元非线性、多元线性，其核心思想都是：多项式拟合；核心方法都是：最小二乘原理。

多元非线性拟合原理

定义：变量 $y$ 与 $m$ 个自变量 $x_j(j=1,2,\cdots,m)$ 都存在线性关系，即构成了多元线性关系：

$y = a_0 + a_1x_1 + a_2x_2 + \cdots + a_mx_m$

每一个自变量 $x_j$ 都有一组(有 $n$ 个)测量值，我们可以记做： $x_{ij}(i=1,2,\cdots,n)$ 。注意： $n>m$ ，因为方程数要比未知数多。变量序号与其观测值序号如下表：

观测值序号\变量序号	$x_1$	$x_2$	$\cdots$	$x_{m-1}$	$x_m$	$y$
1	$x_{11}$	$x_{12}$	$\cdots$	$x_{1(m-1)}$	$x_{1m}$	$y_1$
2	$x_{21}$	$x_{22}$	$\cdots$	$x_{2(m-1)}$	$x_{2m}$	$y_2$
$\vdots$	$\vdots$	$\vdots$	$\vdots$	$\vdots$	$\vdots$	$\vdots$
n	$x_{n1}$	$x_{n2}$	$\cdots$	$x_{n(m-1)}$	$x_{nm}$	$y_n$

如何衡量拟合效果？还是根据最小二乘原理，只不过是未知参数多了一些而已。现在的未知数为 $a_i$ ，即将最小二乘原理运用到包含 $m+1$ 个 $a_i$ 的函数上即可：

$设\quad F = F(a_0,a_1,\cdots,a_m)$

$F = \sum_{i=1}^{n}\left[ y_i - (a_0 + a_1x_{i1} + a_2x_{i2} + \cdots + a_mx_{im} )\right]^2$

对于上式，我们一次对 $F$ 求各个 $a_i$ 的偏导数并令其值为0，可得如下 $m+1$ 个方程：

$\frac{\partial F}{\partial a_0} = 0，\frac{\partial F}{\partial a_1} = 0，\frac{\partial F}{\partial a_2} = 0，\cdots，\frac{\partial F}{\partial a_m} = 0$

最后用于求解的正规方程组为：

$a_0n \quad + a_1\sum_{i=1}^{n}x_{i1} + \quad a_2\sum_{i=1}^{n}x_{i2}\quad + \cdots + a_m\sum_{i=1}^{n}x_{im} \quad\quad = \sum_{i=1}^{n}y_{i}$

$a_0\sum_{i=1}^{n}x_{i1} + a_1\sum_{i=1}^{n}x_{i1}x_{i1} + a_2\sum_{i=1}^{n}x_{i1}x_{i2} + \cdots + a_m\sum_{i=1}^{n}x_{i1}x_{im} = \sum_{i=1}^{n}x_{i1}y_{i}$

$\cdots$

$a_0\sum_{i=1}^{n}x_{im} + a_1\sum_{i=1}^{n}x_{im}x_{i1} + a_2\sum_{i=1}^{n}x_{im}x_{i2} + \cdots + a_m\sum_{i=1}^{n}x_{im}x_{im} = \sum_{i=1}^{n}x_{im}y_{i}$

解上面的线性方程组即可得到一些列系数。

一元非线性拟合

为什么说一元非线性拟合和多元线性拟合关系密切呢？因为为了计算方便，我们会使用变量替换的方法，将一元非线性各阶式子用不同的变量代替，转换成多元非线性，然后用上面的正规方程组求解。

例如：已知有 $n$ 组非线性数据 $(x_i,y_i)(i = 1,2,\cdots,n)$ ，用m阶一元非线性多项式拟合：

$y = a_0 + a_1x + a_2x^2 + \cdots + a_mx^m，n>m$

我们只需要用 $t_j = x^j$ 代换( $j = 1,2,\cdots,m$ )，就可以将一元非线性转换为多元线性：

$y = a_0 + a_1t_1 + a_2t_2 + \cdots + a_mt_m$

用上面多元线性的操作，根据最小二乘原理令关于每个未知系数的偏导数为0，可得正规方程的矩阵的通用形式如下：

图1：一元非线性拟合最终系数的求解方程

关于上面最终求解方法有3个细节要注意：

观测值个数 $n$ 只出现在求和的上限，故不影响矩阵结构；
系数 $a_i$ 的个数，要比多项式最高阶数多1个，即多了前面的 $a_0$ ；
系数矩阵为对称方阵，尺寸比多项式最高阶次 $m$ 多1。

举一个例子：

如何现在有 $n=10$ 个观测数据，想要最高阶数 $m=2$ 的多项式拟合，即：

$y = a_0 + a_1x +a_2x^2$

最终求解矩阵(尺寸： $m+1 = 3$ )就是：

图2：例子的最终求解矩阵

用Matlab实现任意阶非线性拟合

数据自己给，然后输入想要的阶数后，自动求解拟合公式：

% 任意多项式拟合数据点, 当然最好不要超过6阶

clear; clc;

% 这里修改原始观测数据:
xnum = 1:10;
ynum = [3.8 6.3 7.9 8.6 9.2 9.5 9.7 9.9 10.1 10.2];

m = double(input('拟合阶数m = '));
n = length(xnum);

% 等号左边矩阵: 每个元素的幂 = row + col -2
A = zeros(m+1,m+1);
for row = 1:m+1
    for col = 1:m+1
        A(row,col) = sum(xnum.^(row+col-2));
    end
end

% 等号右边矩阵: 一行m+1列
B = zeros(m+1,1);
for num = 1:m+1
    B(num,1) = sum((xnum.^(num-1)).*ynum);
end

% 多项式系数/方程求解:
a = inv(A)*B;

syms x;
ytmp = sym(zeros(1,m+1));
for num = 1:num
    ytmp(num) = x^(num-1);
end
fprintf('拟合结果为:\n');
y = vpa(sum(ytmp.*a'),6)   % 拟合的多项式结果

figure(1);
scatter(xnum,ynum);
hold on;
x = min(xnum):0.1:max(xnum);
y = double(subs(y));
plot(x,y);
grid on;
xlabel('x'); ylabel('y');
title('原点为样本点;实线为多项式拟合结果');

例如5阶效果：

图3：5阶非线性拟合结果

数值分析：多元线性拟合和一元非线性拟合

数值分析：多元线性拟合和一元非线性拟合

前言

多元非线性拟合原理

一元非线性拟合

用Matlab实现任意阶非线性拟合