加权最小二乘估计

前一节中,我们讲解了《最小二乘估计》。在本文中,我们将讲解加权最小二乘估计,这是一种建立在最小二乘估计之上的经典估计方法。

我将全文分为以下三个部分:

  • 绪论
  • 加权最小二乘推导
  • 示例

本文的绝大部分内容从《Optimal State Estimation》一书中翻译而来,并加上了一些自己的学习体会。

绪论

在《最小二乘估计》中,我们假设对所有测量值都有同样的信心。现在假设对其中的一部分测量值,我们对他们的准确度更有信心。这种对所有测量依据信心程度加权的最小二乘估计,就成为加权最小二乘估计。

举例来说,测量一个阻值未知的电阻,我们首先使用一种昂贵的低噪声万用表在清晨精力充沛地测得了一些数据,之后又用一个便宜的万用表在凌晨四点半半死不活地测得了一些数据。显然我们认为前者更加准确一些。后者虽然可能不够精确,但是它们有没有利用价值呢?答案是确定的,我们仍然可以从中获取到一些有用信息。我们永远不要因为数据不准确就把他们抛弃(我们可以加权,星星之火来燎原嘛~)。

加权最小二乘推导

假设x是一个由n个常值元素组成的未知矢量。y是一个由k个元素组成的带有噪声的测量矢量。我们假设y的每个元素都由两部分组成:1. x中每个元素的线性组合。2. 叠加上一个噪声:

left[ {begin{array}{<em>{20}{c}}{{y_1}}\\ vdots \\{{y_k}}end{array}} right] = left[ {begin{array}{</em>{20}{c}}{{H_{11}}}& cdots &{{H_{1n}}}\\ vdots & ddots & vdots \\{{H_{k1}}}& cdots &{{H_{kn}}}end{array}} right]left[ {begin{array}{<em>{20}{c}}{{x_1}}\\ vdots \\{{x_n}}end{array}} right] + left[ {begin{array}{</em>{20}{c}}{{v_1}}\\ vdots \\{{v_k}}end{array}} right]

我们假设每个测量噪声都是零均值且相互独立。测量噪声协方差矩阵为:

begin{array}{l}R = E(v{v^T})\\;;; = left[ {begin{array}{*{20}{c}}{sigma _1^2}& cdots &0\\ vdots & ddots & vdots \\0& cdots &{sigma _k^2}end{array}} right]end{array}

残差的表示与《最小二乘估计》中相同。现在我们构造加权最小二乘估计的代价函数:

J = frac{{varepsilon _{y1}^2}}{{sigma _1^2}} +  cdots  + frac{{varepsilon _{yk}^2}}{{sigma _k^2}}

我们可以将此式与最小二乘估计中的代价函数作对比,后者为:

J = varepsilon _{y1}^2 +  cdots  + varepsilon _{y1}^2

我们发现在加权最小二乘中的代价函数中,残差都多了一个分母,而这个分母为测量噪声的方差,这就是所谓的“加权”。

在进行下一步讨论之前,我们先来仔细看看这个“加权数”:对于测量分量({y_i}),如果它受到剧烈的噪声影响,那么这个噪声的方差就大,因此在加权最小二乘代价函数中这一项的分母就大。带着一个大分母,这一项在代价函数中的重要程度就下降了,也就是说,在下一步对代价函数求极值时,这一项对最终的最优解的影响程度减小了。因此,所谓的加权最小二乘,实际是利用噪声的方差(剧烈程度)进行加权。

下面我们对代价函数进行推导,用矩阵形式展开:

begin{array}{l}J = frac{{varepsilon _{y1}^2}}{{sigma _1^2}} +  cdots  + frac{{varepsilon _{yk}^2}}{{sigma _k^2}}\\;;; = varepsilon _y^T{R^{ - 1}}{varepsilon _y}\\;;; = {(y - Hhat x)^T}{R^{ - 1}}(y - Hhat x)\\;;; = {y^T}{R^{ - 1}}y - {{hat x}^T}{H^T}{R^{ - 1}}y - {y^T}{R^{ - 1}}Hhat x + {{hat x}^T}{H^T}{R^{ - 1}}Hhat xend{array}

现在我们对代价函数J对hat x求偏导,并令其为0来解方程:

begin{array}{l} frac{{partial J}} {{partial hat x}} =  - {y^T}{R^{ - 1}}H + {{hat x}^T}{H^T}{R^{ - 1}}H\\ ;;;;; = 0\\{H^T}{R^{ - 1}}y = {H^T}{R^{ - 1}}Hhat x\\;;;;;;;;;hat x = {({H^T}{R^{ - 1}}H)^{ - 1}}{H^T}{R^{ - 1}}yend{array}

这样,我们就得出了加权最小二乘估计!

这里有一点要注意的是,矩阵R必须是非奇异的,换句话说,测量必须要认为是收到噪声影响的。

在这里,我们直接给出最小二乘估计解作为对比:

hat x = {({H^T}H)^{ - 1}}{H^T}y

示例

让我们回到文章开头,估计阻值的那个例子,我们试图从万用表测量的k个含有噪声的测量中估计一个未知电阻的阻值。在本例中,x是一个标量,我们的测量方程如下:

begin{array}{l}{y_i} = x + {v_i}\\E(v_i^2) = sigma _i^2;;(i = 1, ldots ,k)end{array}

我们将测量写成矩阵形式:

left[ {begin{array}{<em>{20}{c}}{{y_1}}\\ vdots \\{{y_k}}end{array}} right] = left[ {begin{array}{</em>{20}{c}}1\\ vdots \\1end{array}} right]x + left[ {begin{array}{*{20}{c}}{{v_1}}\\ vdots \\{{v_k}}end{array}} right]

测量噪声协方差矩阵为:

R = diag(sigma _1^2, ldots ,sigma _k^2)

直接套用加权最小二乘估计的求解公式:

begin{array}{l}hat x = {({H^T}{R^{ - 1}}H)^{ - 1}}{H^T}{R^{ - 1}}y\\;;; = {left( {[1; cdots ;1]{{left[ {begin{array}{<em>{20}{c}}{sigma _1^2}& cdots &0\\ vdots & ddots & vdots \\0& cdots &{sigma _k^2}end{array}} right]}^{ - 1}}left[ {begin{array}{</em>{20}{c}}1\\ vdots \\1end{array}} right]} right)^{ - 1}} times \\;;;;[1; cdots ;1]{left[ {begin{array}{<em>{20}{c}}{sigma _1^2}& cdots &0\\ vdots & ddots & vdots \\0& cdots &{sigma _k^2}end{array}} right]^{ - 1}}left[ {begin{array}{</em>{20}{c}}{{y_1}}\\ vdots \\{{y_k}}end{array}} right]\\;; = {(sum {frac{1}{{sigma _i^2}}} )^{ - 1}}(frac{{{y_1}}}{{sigma _1^2}} +  cdots  + frac{{{y_k}}}{{sigma _k^2}})end{array}

我们从中可以看书最优解是对测量的加权。每个测量都根据其不确定度的逆(inverse)加权。换句话说,我们凭借自己对测量准确度的知识更加重视部分测量值。如果所有测量的方差相等,那么加权最小二乘估计便退化为最小二乘估计。