针对双变量含有误差的直线拟合问题的方法和设备与流程

本发明一般地涉及测量平差领域。更具体地，本发明涉及用于解决双变量含有误差的直线拟合问题的方法和设备。

背景技术：

直线拟合在各类数据分析或处理中经常遇到，它可以描述为：对于给定的一组观测数据(xi，yi)，(i＝1，2，...，n)，寻找一条最佳的拟合直线y＝kx+b(或直线斜率k和截距b的最佳估计)。对上述问题，如果自变量x的观测值不含误差(仅y的观测值包含误差)，则可直接以y＝kx+b作为观测模型，通过经典的最小二乘平差理论解决。然而，如果x和y的观测值均含有误差，直线拟合问题就变得更加复杂，上述观测模型不再适用。

实际上，双变量存在误差的直线拟合问题，是测量平差中“观测向量、系数矩阵都存在误差”这类问题在直线拟合中的特例。针对这类测量平差问题，golub和vanloan提出了奇异值分解(singularvaluedecomposition，svd)算法，并采用整体最小二乘(totalleastsquares，tls)命名，区别于经典的最小二乘估计，tls适用于观测数据等精度、独立的情况；针对实际问题中观测数据常常不等精度、相关的情况，schaffrin和wieser提出了加权整体最小二乘(weightedtotalleastsquares，wtls)方法，许多研究者也相继提出了一些迭代算法。

虽然双变量有误差情况下的直线拟合问题可以通过wtls的理论和模型正确解决，但其数学模型和迭代格式复杂，也没有充分利用这一特例问题的特点和经典平差理论的优势。

技术实现要素：

为至少解决上述技术问题，在一个方面中，本发明提供了一种用于解决双变量含有误差的直线拟合问题的方法和设备，其中所述方法包括：

建立间接平差的函数模型及随机模型；

初始化参数x⁰并求解参数改正量

更新参数估计值并迭代求解直至收敛，其中所述参数估计值为以及

执行精度评定。

在一个实施例中，所述建立间接平差的函数模型及随机模型的过程包括：取直线的斜率k、截距d及n个观测点的x坐标为参数，即取n个观测点的x坐标、y坐标为观测量，即根据间接平差模型，将观测量对参数在处线性化，从而得出间接平差的误差方程中a和l的表达式如下：

误差方程中为待求解的参数改正量，取单位权方差为1，则上述观测量的权阵为：

在另一个实施例中，所述初始化参数x⁰并求解参数改正量的过程包括：首先利用(x1，y1)、(xn，yn)求解k⁰、b⁰如下式(3)：

并取完成x⁰初始化。然后根据式(1)计算a和l，由下式(4)计算参数改正量

在另一个实施例中，所述更新参数估计值并迭代求解直至收敛的过程包括：判断参数改正量(向量)的第1、2元素，即直线斜率改正量截距改正量的绝对值与收敛门限tk、tb的大小：

1、若且退出迭代；

2、否则，取作为新的初始化值x⁰，根据式(1)计算a和l，根据式(4)计算对应新的初始化值的参数改正量当迭代结束时，得到参数的估计值：

该参数(向量)的第1、2元素，即直线斜率、截距的估计值。

在另一个实施例中，所述执行精度评定的过程包括：

根据间接平差理论和公式，按照下式(6)计算拟合后单位权方差：

按照下式(7)计算参数的中误差：

∑＝σ0(a^tpa)^-1(7)

该矩阵对角线的第1、2元素，即直线斜率、截距估计值的中误差。

在另一个方面中，本发明还提供了一种用于解决双变量含有误差的直线拟合问题的设备，所述设备包括：

建立模块，其配置成用于建立间接平差的函数模型及随机模型；

初始化模块，其配置成用于初始化变量x⁰并求解参数改正量

更新求解模块，其配置成用于更新参数估计值并迭代求解直至收敛，其中所述参数估计值为以及

精度评定模块，其配置成用于执行精度评定。

在一个实施例中，所述建立模块还配置成用于执行：

取直线的斜率k、截距d及n个观测点的x坐标为参数；

取n个观测点的x坐标、y坐标为观测量

根据间接平差模型，将观测量对参数在处线性化，从而得出间接平差的误差方程中a和l的表达式如下：

误差方程中为待求解的参数改正量，取单位权方差为1，则上述观测量的权阵为：

在另一个实施例中，所述初始化模块还配置成用于执行：

首先利用(x1，y1)、(xn，yn)求解k⁰、b⁰如下式(3)：

并取完成x⁰初始化。

然后，根据式(1)计算a和l，由下式(4)计算参数改正量

在另一个实施例中，所述更新求解模块还配置成用于执行：

判断参数改正量(向量)的第1、2元素，即直线斜率改正量截距改正量的绝对值与收敛门限tk、tb的大小：

1、若且退出迭代；

2、否则，取作为新的初始化值x⁰，根据式(1)计算a和l，根据式(4)计算对应新的初始化值的参数改正量当迭代结束时，得到参数的估计值：

该参数(向量)的第1、2元素，即直线斜率、截距的估计值。

在另一个实施例中，所述执行模块还配置成用于执行：

根据间接平差理论和公式，按照下式(6)计算拟合后单位权方差：

按照下式(7)计算参数的中误差：

∑＝σ0(a^tpa)^-1(7)

该矩阵对角线的第1、2元素，即直线斜率、截距估计值的中误差。

本发明基于经典最小二乘平差理论，在考虑双变量含有误差的情况下，提出了新的直线拟合观测模型，通过重新构建直线拟合的函数模型，解决了双变量含有误差情况下直线斜率和截距的最佳估计问题，并完成精度评定。

该方法采用间接平差方法，利用直线拟合问题的特殊性及经典平差理论的优势，获得直线斜率和截距的最佳估计，并完成精度评定，避免了wtls方法的复杂模型和迭代格式，模型清晰，数学形式简单，易于编程计算。该方法可用于解决测绘、物理、航天等诸多领域的双变量含有误差的直线拟合问题，量化分析双变量之间的线性关系，解决相应领域的具体实际技术问题。

附图说明

通过阅读仅作为示例提供并且参考附图进行的以下描述，将更好地理解本发明及其优点，其中：

图1是示出根据本发明的用于解决双变量含有误差的直线拟合问题的方法的流程图；

图2是示出根据本发明的用于解决双变量含有误差的直线拟合问题的设备组成框图；以及

图3是示出采用本发明的方法对一组观测数据进行直线拟合的结果。

具体实施方式

本发明针对双变量含有误差的直线拟合问题，提出一种新的直线拟合方法，该方法采用间接平差原理，获得了直线斜率和截距的最佳估计，并完成精度评定。

下面结合附图和具体实施方式对本发明进行详细说明。

图1是根据本发明的用于解决双变量含有误差的直线拟合问题的方法的流程框图。如图1所示，本发明所述方法100首先开始执行步骤101，在该步骤中，建立间接平差的函数模型及随机模型；接着方法100进入到102步骤，在该步骤中，初始化参数x⁰并求解参数改正量然后方法100前进到103步骤，在该步骤中，更新参数估计值并迭代求解直至收敛，其中所述参数估计值为最后，方法100执行精度评定104步骤。

假设一组观测数据为(xi，yi)(i＝1，2，...，n)，这些观测数据的方差记为

在一个实施例中，建立间接平差的函数模型及随机模型的过程包括：取直线的斜率k、截距d及n个观测点的x坐标为参数，即取n个观测点的x坐标、y坐标为观测量，即

根据间接平差模型，将观测量对参数在处线性化，从而间接平差的误差方程中a和l的表达式如下：

误差方程中为待求解的参数改正量，取单位权方差为1，则上述观测量的权阵为：

在另一个实施例中，初始化x⁰并求解参数改正量102的过程包括：初始化x⁰，并求解参数改正量利用(x1，y1)、(xn，yn)求解k⁰、b⁰如下：

并取完成x⁰初始化。

根据式(1)计算a和l，由下式(4)计算参数改正量

在另一个实施例中，更新参数估计值并迭代求解直至收敛103的过程包括：判断参数改正量(向量)的第1、2元素，即直线斜率改正量截距改正量的绝对值与收敛门限tk、tb的大小：

1、若且退出迭代；

2、否则，取作为新的初始化值x⁰，根据式(1)计算a和l，根据式(4)计算对应新的初始化值的参数改正量

当迭代结束时，得到参数的估计值：

该参数(向量)的第1、2元素，即直线斜率、截距的估计值。

在另一个实施例中，执行精度评定104的过程包括：

根据间接平差理论和公式，按照下式(6)计算拟合后单位权方差：

按照下式(7)计算参数的中误差：

∑＝σ0(a^tpa)^-1(7)

该矩阵对角线的第1、2元素，即直线斜率、截距估计值的中误差。

图2是示出根据本发明的用于解决双变量含有误差的直线拟合问题的设备组成框图。如图2所示，本发明所述设备200包括：建立模块201，其配置成用于建立间接平差的函数模型及随机模型；初始化模块202，其配置成用于初始化变量x⁰并求解参数改正量更新求解模块203，其配置成用于更新参数估计值并迭代求解直至收敛，其中所述参数估计值为以及精度评定模块204，其配置成用于执行精度评定。

在一个实施例中，所述建立模块201还配置成用于执行：

取直线的斜率k、截距d及n个观测点的x坐标为参数

取n个观测点的x坐标、y坐标为观测量

根据间接平差模型，将观测量对参数在处线性化，从而得出间接平差的误差方程中a和l的表达式如下：

误差方程中为待求解的参数改正量，取单位权方差为1，则上述观测量的权阵为：

在另一个实施例中，所述初始化模块202还配置成用于执行：

首先利用(x1，y1)、(xn，yn)求解k⁰、b⁰如下式(3)：

并取完成x⁰初始化。

然后，根据式(1)计算a和l，由下式(4)计算参数改正量

在另一个实施例中，所述更新求解模块203还配置成用于执行：

判断参数改正量(向量)的第1、2元素，即直线斜率改正量截距改正量的绝对值与收敛门限tk、tb的大小：

1、若且退出迭代；

2、否则，取作为新的初始化值x⁰，根据式(1)计算a和l，根据式(4)计算对应新的初始化值的参数改正量当迭代结束时，得到参数的估计值：

该参数(向量)的第1、2元素，即直线斜率、截距的估计值。

在另一个实施例中，所述执行模块204还配置成用于执行：

根据间接平差理论和公式，按照下式(6)计算拟合后单位权方差：

按照下式(7)计算参数的中误差：

∑＝σ0(a^tpa)^-1(7)

该矩阵对角线的第1、2元素，即直线斜率、截距估计值的中误差。

下面结合图3的结果示意图具体说明本发明所述方法的实现原理。为方便与其他方法比较，采用与文献schaffrinbandwiesera.onweightedtotalleast-squaresadjustmentforlinearregression[j].journalofgeodesy，2008，82(7)：415-421.中相同的数据，如下：

基于以上数据，按照如下流程进行直线拟合。

步骤1：建立间接平差的函数模型及随机模型101。根据以上数据，n＝10，即将观测量对参数在处线性化，从而a和l如下：

上述观测量的权阵为：

上式中，diag()表示以括号中的元素为对角线上的元素构造对角阵。

步骤2：初始化x⁰并求解参数改正量102。利用(x1，y1)、(x10，y10)求解k⁰、b⁰如下：

并取完成x⁰初始化。

根据式(9)计算a和l，并计算参数改正量

步骤3：更新参数估计值并迭代求解直至收敛103。设收敛门限tk、tb均为1×10^-10，由于0.130852989、0.503485017均大于1×10^-10，需取作为新的初始化值x⁰，即

x⁰＝[-0.4637416065.396514983…]^t(19)

根据式(9)计算a和l，同时计算参数改正量并判断收敛情况。如此迭代往复，直至得到第13次计算所得的参数改正量

此时满足收敛条件，退出迭代，进而

由此得出直线斜率、截距的估计值分别为-0.480533407和5.479910224。

步骤4：执行精度评定，其中包括计算拟合后单位权方差：

按照下式计算参数的中误差：

由此直线斜率、截距估计值的中误差分别为0.07、0.36。

根据本发明计算所得结果与参考文献、真值的对比如下：

图3是采用本发明的方法对如上一组观测数据进行直线拟合的结果与真值进行比较的示意图。从以上对比及图3可以看出，本发明与参考文献所得结果完全一致，且逼近真值，说明本发明的模型及方法正确。

本发明的针对双变量含有误差的直线拟合问题的方法可以通过在计算机可读的记录介质以计算机可读代码来体现。计算机可读记录介质包括存储可通过计算机系统解读的数据的所有种类的记录介质。该记录介质例如可以包括但不限于只读存储器(rom，“readonlymemory”)、随机存取存储器(ram，“randomaccessmemory”)、磁盘、磁盘、光盘、闪存等。进一步，这些计算机可读的记录介质可以通过通信网络(包括计算机通信网络、蜂窝通信网络或局域域通信网络)在各个通信实体之间传播或扩散，从而也可以通过任意的方式来运行存储在计算机可读存储介质上的计算机可读指令或计算机可执行代码。

虽然本发明所实施的方式如上，但所述内容只是为便于理解本发明而采用的实施例，并非用以限定本发明的范围和应用场景。任何本发明所述技术领域内的技术人员，在不脱离本发明所揭露的精神和范围的前提下，可以在实施的形式上及细节上作任何的修改与变化，但本发明的专利保护范围，仍须以所附的权利要求书所界定的范围为准。