MATLAB在数据误差处理中的应用

一 前言
 随着国民经济的迅速发展，大量的数据需要处理，误差理论和数据处理的任务也越来越重，传统的手算以及传统的计算器等工具已不能满足需要。另一方面，计算机在我们的日常生活中却日益普及，显然，运用计算机进行数据处理已是大势所趋。
 MATLAB是美国MathWorks公司推出的一种简洁方便的工程计算语言，自从问世就以其友好的用户界面和多种功能深受各方面用户的欢迎。
 测量数据的数据处理和数据分析涉及到最小二乘法、回归分析、曲线拟合以及线性方程组的求解等内容，而这些正是MATLAB的强项。另外，通过MATLAB强大的图形功能，我们还能方便地将数据图形化，从而进行直观地分析处理数据。
 二 几个基本概念
 1、误差：
 若令测量误差为，测得值为x,真值为，则有
 或  （2-1）
 由于实际应用中真值一般是无从知道或无法确定的，所以，在统计学中，常以测量次数足够大时的测得值的算术平均值近似代替真值。
 2、算术平均值：
 对一真值为的物理量进行等精度的n次测量，得n个测得值，它们都含有随机误差，统称真差。我们常以算术平均值作为n次测量的结果，即
   （2-2）
 3、残差v：
 各测得值对其算术平均值的误差量叫做残余误差，简称残差，即
   （2-3）
 4、标准差（标准偏差）：
 在计量学中，常用标准差来评定测得值的精度，即
   （2-4）
 式中：：真差（随机误差）；
 n：测量次数。
 但在实际应用中，真差往往是不可知的，而常根据有限个测量值的残差v来求取随机测量误差方差的估计值，开方，得
   （2-5）
 式2-5称为贝塞尔（Bessel）公式，称为试验标准差，即是标准差的估计值。
 5、随机误差的正态分布：
 正态分布是随机误差的一种重要分布。实践表明，在大多数情况下，在测量过程中产生的误差服从正态分布。
 正态分布的分布曲线如图1所示，
 
 其分布密度函数为
   （2-6）
 式中，y：概率密度；
 x：随机变量；
 ：标准差；
 ：理论均值或随机变量x的数学期望。
 因被测量的真值无法知道，对连续型随机函数，可将理论均值看作真值，故式2-6可写作
   （2-7）
 若用代替，则分布密度函数又可化为
  （2-8）
 式2-8说明，测量次数足够大时，正态分布方程式同样适用于残差v。
 6、非等精度测量的加权平均值及其精度参数：
 “权”即各组测量结果相对的可信赖程度，一般用符号p代表“权”，所以求取加权平均值可使用下式
   （2-9）
 而各组测量的“权”，与各组测量结果的方差成反比，即
   （2-10）
 单位权化以后所得的单位权的标准差为
   （2-11）
 式中：
 m：测量组数。
 而加权平均值的标准差为
   （2-12）
 三 在计算几个基本的数字特征中的应用
 1、求算术平均值：
 计算一组数据的算术平均值，使用mean函数，其语法格式为：
 m=mean(x)
 x为所求的一组数据组成的行向量。
测量一个长度10次，所得结果如表1，求数据的算术平均值：
表1
序号 1 2 3 4 5 6 7 8 9 10
长度(mm) 25.125 25.126 25.127 25.128 25.129 25.130 25.131 25.132 25.133 25.134
 程序如下，可得结果为＝25.1295。
 >> y=25.125:.001:25.134;
 >> m=mean(y)
 m =
    25.1295
 2、求残差v：
 计算一组数据样本的程序十分简单，故MATLAB中没有相应的子程序供调用，但我们可以用下面的程序进行求解（设m是数据样本的算术平均值）：
 例2．求例1中的数据样本的残差：
 程序如下：
 >> y=25.125:.001:25.134;
 >> m=mean(y);
 >> vi=y-m
 vi =
    -0.0045   -0.0035   -0.0025   -0.0015   -0.0005    0.0005    0.0015    0.0025    0.0035    0.0045
 所得vi即为所求的残差。
 3、求标准差：
 如二.4所述，计算一组数据的标准差，常用计算试验标准差代替，此时使用std函数，其格式为
 ＝std(x)
 x为数据样本组成的一组行向量。
 例3．计算例1中的数据的标准差：
 程序如下，可得结果为＝0.0030。
 >> y=25.125:.001:25.134;
 >> s=std(y)
 s =
     0.0030
 4、正态分布的随机误差的一些参数的求法：
 数据样本的随机误差多服从正态分布，用normstat函数求正态分布的均值和方差，其语法格式为：
 [m,v]=normstat(mu,sigma)
 5、计算非等精度测量的加权平均值及其精度参数：
 例4．1m米尺由3位观测者测量，其结果如表2，求加权平均值及标准差：
 表2
组别 一 二 三
（mm） 1000.045 1000.015 1000.060
（m） 5 20 10
 程序如下：
 >> format long
 >> sig=[5 20 10];
 >> mx=[1000.045 1000.015 1000.060];
 >> p=1./sig.^2
 p =
     0.0400    0.0025    0.0100
 >> p=p.*400
 p =
     16     1     4
 >> xp=sum(mx.*p)/sum(p)  %求加权平均值
 xp =
     1.000046428571429e+003 
 >> sis=sqrt(sum(p.*(mx-xp).^2)/(3-1))  %求单位权组的标准差
 sis =
    0.02964070560175
 >> simx=sis/(sqrt(sum(p)))  %求加权平均值的标准差
 simx =
    0.00646813224152
 四 在使用最小二乘法时的应用
 众所周知，最小二乘法在数据处理中具有无法取代的重要地位。最小二乘法既可处理满足线性函数关系的数据样本，也可以处理满足非线性函数关系的数据样本。
 1、线性函数的最小二乘法处理：
 已知数据样本符合线性函数关系，即：y=ax+b，测得的数据样本为长度相等的x,y向量。在MATLAB中通常使用矩阵除法来求解：
 设矩阵A、c、y如下：
 、、
 则问题可化为解线性方程：cA=y，在MATLAB中可用c=A\y进行求解，求得列向量c即可得出系数a=c(1,1)、b=(2,1)，然后得出线性函数关系。
 例5．为研究20mm轴的几何形状误差，在40mm长度内选5个断面测得直径偏差如表3，试确定沿长度方向形状误差的规律。
 表3 数据表
被测断面距端面距离/mm 2 10 20 30 40
直径偏差/m +3 +5 +8 +15 +18
 解：先将（，）图示于图2上（图中黑色“＋”所示），
 经初步分析即知误差呈线性规律。
 设此规律的线性方程为：
 然后在MATLAB中用最小二乘法线性拟合可得近似y、x值为：y=0.4185； x=1.2617。
 程序如下：
 >> li=[2 10 20 30 40];
 >> dd=[3 5 8 15 18];
 >> a=[li;ones(size(li))]';
 >> y=dd';
 >> c=A\y
 c =
     0.4185
     1.2617
 >> A=c(1,1),b=c(2,1)
 A =
     0.4185
 b =
     1.2617
 所以，所求的规律近似为：
 图示于图2（图中黑色实线）。
 尽管MATLAB中没有直接供调用的最小二乘法处理系统函数，但我们可以自己直接编写.m文件来供调用，文件如下：
 function [a,b]=l