考虑一组
(
x
1
,
y
1
)
,
(
x
2
,
y
2
)
,
…
,
(
x
m
,
y
m
)
{\displaystyle (x_{1},y_{1}),(x_{2},y_{2}),\dots ,(x_{m},y_{m})}
共
m
{\displaystyle m}
个数据点以及曲线(模型函数)
y
^
=
f
(
x
,
β
)
{\displaystyle {\hat {y}}=f(x,{\boldsymbol {\beta }})}
。该曲线同时取决于x 与
β
=
(
β
1
,
β
2
,
…
,
β
n
)
{\displaystyle {\boldsymbol {\beta }}=(\beta _{1},\beta _{2},\dots ,\beta _{n})}
共n 个参数(满足
m
≥
n
{\displaystyle m\geq n}
)。目标是找到在最小二乘意义上与数据点拟合最好的曲线所对应的参数
β
{\displaystyle {\boldsymbol {\beta }}}
,即最小化平方和
S
=
∑
i
=
1
m
r
i
2
,
{\displaystyle S=\sum _{i=1}^{m}r_{i}^{2},}
其中残差 ri 的定义为
r
i
=
y
i
−
f
(
x
i
,
β
)
,
(
i
=
1
,
2
,
…
,
m
)
.
{\displaystyle r_{i}=y_{i}-f(x_{i},{\boldsymbol {\beta }}),\qquad (i=1,2,\dots ,m).}
S 取最小值 时的梯度 为零。由于模型包含n 个参数,因此可得到n 个梯度方程:
∂
S
∂
β
j
=
2
∑
i
r
i
∂
r
i
∂
β
j
=
0
(
j
=
1
,
…
,
n
)
.
{\displaystyle {\frac {\partial S}{\partial \beta _{j}}}=2\sum _{i}r_{i}{\frac {\partial r_{i}}{\partial \beta _{j}}}=0\quad (j=1,\ldots ,n).}
在非线性系统中,偏导数
∂
r
i
∂
β
j
{\textstyle {\frac {\partial r_{i}}{\partial \beta _{j}}}}
同时是自变量x 和参数
β
{\displaystyle {\boldsymbol {\beta }}}
的函数,因此这些梯度方程通常没有封闭解。因而必须为参数选择初始值用以迭代求解。迭代表达式为
β
j
≈
β
j
k
+
1
=
β
j
k
+
Δ
β
j
.
{\displaystyle \beta _{j}\approx \beta _{j}^{k+1}=\beta _{j}^{k}+\Delta \beta _{j}.}
其中,k 是迭代次数,
Δ
β
{\displaystyle \Delta {\boldsymbol {\beta }}}
则是偏移向量。每次迭代时,使用关于
β
k
{\displaystyle {\boldsymbol {\beta }}^{k}}
的一阶泰勒级数 展开以线性化模型:
f
(
x
i
,
β
)
≈
f
(
x
i
,
β
k
)
+
∑
j
∂
f
(
x
i
,
β
k
)
∂
β
j
(
β
j
−
β
j
k
)
=
f
(
x
i
,
β
k
)
+
∑
j
J
i
j
Δ
β
j
.
{\displaystyle f(x_{i},{\boldsymbol {\beta }})\approx f(x_{i},{\boldsymbol {\beta }}^{k})+\sum _{j}{\frac {\partial f(x_{i},{\boldsymbol {\beta }}^{k})}{\partial \beta _{j}}}\left(\beta _{j}-\beta _{j}^{k}\right)=f(x_{i},{\boldsymbol {\beta }}^{k})+\sum _{j}J_{ij}\,\Delta \beta _{j}.}
雅可比矩阵 J 是常数、自变量与参数的函数,因此每次迭代时的J 并不固定。对线性化模型而言,
∂
r
i
∂
β
j
=
−
J
i
j
,
{\displaystyle {\frac {\partial r_{i}}{\partial \beta _{j}}}=-J_{ij},}
残差的表达式则为
Δ
y
i
=
y
i
−
f
(
x
i
,
β
k
)
,
{\displaystyle \Delta y_{i}=y_{i}-f(x_{i},{\boldsymbol {\beta }}^{k}),}
r
i
=
y
i
−
f
(
x
i
,
β
)
=
(
y
i
−
f
(
x
i
,
β
k
)
)
+
(
f
(
x
i
,
β
k
)
−
f
(
x
i
,
β
)
)
≈
Δ
y
i
−
∑
s
=
1
n
J
i
s
Δ
β
s
.
{\displaystyle r_{i}=y_{i}-f(x_{i},{\boldsymbol {\beta }})=\left(y_{i}-f(x_{i},{\boldsymbol {\beta }}^{k})\right)+\left(f(x_{i},{\boldsymbol {\beta }}^{k})-f(x_{i},{\boldsymbol {\beta }})\right)\approx \Delta y_{i}-\sum _{s=1}^{n}J_{is}\Delta \beta _{s}.}
将上述表达式代入梯度方程,可以得到
−
2
∑
i
=
1
m
J
i
j
(
Δ
y
i
−
∑
s
=
1
n
J
i
s
Δ
β
s
)
=
0
,
{\displaystyle -2\sum _{i=1}^{m}J_{ij}\left(\Delta y_{i}-\sum _{s=1}^{n}J_{is}\ \Delta \beta _{s}\right)=0,}
以上方程可化简为n 个联立的线性方程,称为正规方程 (normal equations):
∑
i
=
1
m
∑
s
=
1
n
J
i
j
J
i
s
Δ
β
s
=
∑
i
=
1
m
J
i
j
Δ
y
i
(
j
=
1
,
…
,
n
)
.
{\displaystyle \sum _{i=1}^{m}\sum _{s=1}^{n}J_{ij}J_{is}\ \Delta \beta _{s}=\sum _{i=1}^{m}J_{ij}\ \Delta y_{i}\qquad (j=1,\dots ,n).}
正规方程可用矩阵表示法写成
(
J
T
J
)
Δ
β
=
J
T
Δ
y
.
{\displaystyle \left(\mathbf {J} ^{\mathsf {T}}\mathbf {J} \right)\Delta {\boldsymbol {\beta }}=\mathbf {J} ^{\mathsf {T}}\ \Delta \mathbf {y} .}
上述方程是使用高斯-牛顿算法 求解非线性最小二乘问题的的基础。
需要注意的是雅可比矩阵定义中导数的符号约定。某些文献中的J 可能与此处的定义相差一个负号。
权重扩展
不同数据点(观测结果)的可靠性并不一定相同,此时可使用加权平方和
S
=
∑
i
=
1
m
W
i
i
r
i
2
.
{\displaystyle S=\sum _{i=1}^{m}W_{ii}r_{i}^{2}.}
权重矩阵W 是一个对角矩阵 ,理想情况下每个权重系数应等于观测误差方差 的倒数。[ 1] 此时,正规方程可扩展为
(
J
T
W
J
)
Δ
β
=
J
T
W
Δ
y
.
{\displaystyle \left(\mathbf {J} ^{\mathsf {T}}\mathbf {WJ} \right)\Delta {\boldsymbol {\beta }}=\mathbf {J} ^{\mathsf {T}}\mathbf {W} \ \Delta \mathbf {y} .}