MATLAB数据误差应用

1、MATLAB在数据误差处理中的应用（于海艳 12010245217 2010级通信班）【摘要】伴随着越来越多的数据的产生，数据处理的任务越来越重，本文就计算机软件MATLAB在数据处理的应用作一大只介绍。文章首先说明了误差理论与数据处理中的几个基本概念，然后详细介绍了数据样本误差请定的几个数字特征在MATLAB中计算方法，接着介绍在的数据处理中广为应用的最小二乘法以及其在MATLAB中的计算过程，文章的最后又大体说明了MATLAB在数据样本的回归分析与利用经验公式求解数据规律中的应用。【关键字】MATLAB 误差 数据处理 正态分布 最小二乘法 回归分析 经验公式一、问题的提出随着国民经济的迅

2、速发展，大量的数据需要处理，误差理论和数据处理的任务也原来越重，传统的手算以及计算器等工具已不能满足需要。另一方面，计算机在我们的日常生活中的越来越普及，显然，运用计算机惊醒数据处理意识大势所趋。 MATLAB是美国MathWorks公司推出的一种简洁方便的工程计算语言，自从其问世就以其友好的用户界面和多种功能深受各方面欢迎。 测量数据的数据处理和数据分析涉及到最小二乘法、回归分析、曲线拟合以及线性方程组的求解内容，而这些正是MATLAB的强项，另外，通过MATLAB强大的图形功能，我们还能方便地将数据图形化，从而进行直观地分析处理数据。二、几个基本概念1、 误差在测量中，误差表示测得值与真值

3、之差，若令测量误差为，测得值为x,真值为，则有=x-x0或x0=x- (1)由于实际应用中真值一般是无从知道或无法确定的，所以，在统计学中，常以测量次数足够大时的测得值的算术平均值近似代替真值。2、 算术平均值 对一真值为x0的物理量进行等精度的n次测量，得n个测得值x1,x2,L,xn，它们都含有随机误差1,2,n，统称真差。我们常以算术平均值作为n次测量的结果，即=x1+ x2+L+ xn)= （2）3、残差v  各测得值对其算术平均值的误差量叫做残余误差，简称残差，即   v= - (3)4、标准差（标准偏差）  在计量学中，常用标准差来评定测得值

4、的精度，即   = (n) （4）  式中：真差（随机误差）；n：测量次数。但在实际应用中，真差往往是不可知的，而常根据有限个测量值的残差v来求取随机测量误差方差的估计值，开方，得   （5） 式2-5称为贝塞尔（Bessel）公式，称为试验标准差，即是标准差的估计值。 5、随机误差的正态分布：  正态分布是随机误差的一种重要分布。实践表明，在大多数情况下，在测量过程中 产生的误差服从正态分布。 图一的程序如下>> x=0:0.02:5;y=1/(.5*sqrt(2*pi)*exp(x*2.5).2/(

5、2*.52); plot(x,y) ylim0,1;xlim0,5 xlabel('x'),ylabel('y')正态分布的分布曲线如图1所示， 其分布密度函数为   y=f(x)= （6） 式中，y：概率密度； x：随机变量； ：标准差； u：理论均值或随机变量x的数学期望。 因被测量的真值无法知道，对连续型随机函数，可将理论均值看作真值，故式2-6可写作   y=f()= （7） 若用代替u，则分布密度函数又可化为  y=f(v)= （8）

6、式2-8说明，测量次数足够大时，正态分布方程式同样适用于残差v。 6、非等精度测量的加权平均值及其精度参数：  “权”对的可信赖程度，一般用符号p代表“权”，所以求取加权平均值可使用下式 （9） 而各组测量的“权”，与各组测量结果的方差成反比，即   P （10） 单位权化以后所得的单位权的标准差为   = （11） m：测量组数。 而加权平均值的标准差为   （12） 三 在计算几个基本的数字特征中的应用 1、求算术平均值 计算一组数据的算术平均值，使用

7、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:0.001:25.134;%将初值为25.125，步长是0.001，末值是25.134的行向量赋值给ym=mean(y)%运用mean

8、函数求y的算术平均值m = 25.1295 %得到m的值 2、求残差v  计算一组数据样本的程序十分简单，故MATLAB中没有相应的子程序供调用，但我们可以用下面的程序进行求解（设m是数据样本的算术平均值）：例2求例1中的数据样本的残差：表一：计算一组数据样本序号12345678910长度（mm）25.12525.12625.12725.12825.12925.13025.13125.13225.13325.134 程序如下： >> y=25.125:0.001:25.134;%将初值是25,125，步长是0.001，末值是25.134的行向

9、量赋值给ym=mean(y);%运用mean函数求y的算术平均值vi=y-m %求vi的表达式vi = Columns 1 through 9 -0.0045 -0.0035 -0.0025 -0.0015 -0.0005 0.0005 0.0015 0.0025 0.0035 Column 10 0.0045>> y=25.125:.001:25.134;% 所得vi即为所求的残差。3、求标准差  计算一组数据的标准差，常用计算试验标准差代替，此时使用std函数，其格式为  std(x) x为数据样本组成的一组行向量。 例3计算例1中的数据

10、的标准差：  程序如下，可得结果为s0.0030。>> y=25.125:.001:25.134; %将初值是25.125，步长是0.001，末值是25.134赋值给y s=std(y) %计算数据的标准差ss =0.0030 %得到标准差的值4、正态分布的随机误差的一些参数的求法：  数据样本的随机误差多服从正态分布，用normstat函数求正态分布的均值和方差，其语法格式为：  m,v=normstat(mu,sigma) 5、计算非等精度测量的加权平均值及其精度参数：  例41m米尺由3位观测者测量，其结果如表2，求加权平均值

11、及标准差：表二：测量数据及标准差组别 一 二 三（mm） 1000.045 1000.015 1000.060 5 20 10 >> format long sig=5 20 10;%定义一个矩阵 mx=1000.045 1000.015 1000.060;% 将平均值以矩阵形式赋值给mx p=1./sig.2 p=p.*400xp=sum(mx.*p)/sum(p)% 求加权平均值 sis=sqrt(sum(p.*(mx-xp).2)/(3-1)%求单位权组的标准差 simx=sis/(sqrt(sum(p)%求加权平均值的标准差p = 0.0400 0.00

12、25 0.0100p = 16 1 4xp = %得到加权平均值 1.0000e+003sis = 0.0296simx = %得到加权平均值的标准差 0.0065 四 、在使用最小二乘法时的应用  众所周知，最小二乘法在数据处理中具有无法取代的重要地位。最小二乘法既可处理满足线性函数关系的数据样本，也可以处理满足非线性函数关系的数据样本。 1、线性函数的最小二乘法处理：  已知数据样本符合线性函数关系，即：y=ax+b，测得的数据样本为长度相等的x,y向量。在MATLAB中通常使用矩阵除法来求解：设矩阵A、c、y如下：  A=、c=、y=则问

13、题可化为解线性方程：cA=y，在MATLAB中可用c=Ay进行求解，求得列向量c即可得出系数a=c(1,1)、b=(2,1)，然后得出线性函数关系。例5为研究20mm轴的几何形状误差，在40mm长度内选5个断面测得直径偏差如表3，试确定沿长度方向形状误差的规律 表三：被测断面的距离及偏差被测断面距面距离/mm210203040直径偏差+3+5+8 +15 +18图二 例5 数据图示 经初步分析即知误差呈线性规律。 设此规律的线性方程为： 然后在MATLAB中用最小二乘法线性拟合可得近似y、x值为：y=0.4185； x=1.2617。 程序如下:li=2

14、 10 20 30 40;%给长度赋值dd=3 5 8 15 18;plot(li,dd,'r+') % 利用plot函数画出函数图象hold on;%保持原图像ylim(gca),0 42）xlim(gca,0 20)% 具有两个元素的数值向量p=polyfit(li,dd,1); %得到最小二乘拟合多项式的系数xx=0:0.2:42; %初值是0，步长为0.2，末值是42的行向量plot(xx,y) 所以，所求的规律近似为： 图示于图2（图中黑色实线）。 尽管MATLAB中没有直接供调用的最小二乘法处理系统函数，但我们可以自己直接编写.m文件来

15、供调用.文件保存为lsline.m，即可供调用。调用程序如下：li=2 10 20 30 40; dd=3 5 8 15 18; a,b=lsline(li,dd) 利用函数文件求系数 a = 0.4185 b = 1.2617 %得到两个系数 即可得出所求的线性规律。 2、非线性函数的最小二乘法处理：  MATLAB中非线性最小二乘的处理使用nlinfit函数，下面我们通过一个例子来介绍它的使用方法。例6.在化工生产中获得的氯气的等级y随生产时间x下降，已知在x>=8时，y与x之间有如下的非线性模型：现收集了10组数据，如表4：表四：y随x非线性变化 X

16、8 10 12 14 16 18 20 22 24 26 y0.490.480.460.430.430.450.410.400.40 0.40x 8 10 12 14 16 18 20 22 24 26 y 0.49 0.48 0.46 0.43 0.43 0.45 0.41 0.40 0.40 0.40  要求利用该数据样本求a、b的值，以确定模型。  首先定义非线性函数的.m文件model.mFunction yy=model(be,x) %定义函数文件modela=be(1);b=be(2);%两个系数yy=a+(0.49-a)*exp(-b*(x-8); %函数表达

17、式值 然后在命令窗口中输入程序：>> x=8 10 12 14 16 18 20 22 24 26; >> y=0.49 0.48 0.46 0.43 0.43 0.45 0.41 0.40 0.40 0.40; %¸定义氯气等级y时间x >> be=0 0; >> befit=nlinfit(x,y,'model',be)%利用nlinfit对函数进行最小二乘法处理 befit = 0.3584 0.0692%得到处理后结果 五 、在回归分析与经验公式中的应用  在日常生活中，人们常应用试

18、验的方法，寻找出数据样本之间的相互关系。但是通常使用的方法往往不能深刻反映变量间的内在关系，而应用经验公式却能充分表达数据样本各变量之间的变化规律，而且便于从理论上作进一步的研究。回归分析法就是应用数理统计的方法，对数据样本进行分析和处理，从而得出反映各变量间相互关系的经验公式，这就是回归方程。 1、一元线性经验公式  一元线性经验公式是指自变量x与因变量y存在的线性变化的规律，其形式为y=ax+b式中：a、b即为需要由数据样本确定的回归参数。在MATLAB中，我们通常使用一次曲线拟合的方法来求解回归参数，曲线拟合的命令语法格式为p=polyfit(x,y,1) %得到最小

19、二乘拟合多项式的系数 求得的p为向量a,b。例7. 用X光机检查镁合金铸件内部缺陷时，为获得最佳灵敏度，透视电压y应随被透视件的厚度x而改变，试验中获得了如表5所示的数据表，试用其中数据确定其一元线性经验公式： 表五：y随x的变化关系x/mm 13 13 14 15 16 18 20 22 24 26y/mm52 55 58 61 65 70 75 80 85 91x/mm 12 13 14 15 16 18 20 22 24 26 y/kV 52 55 58 61 65 70 75 80 85 91 程序如下，得所求公式为：y=2.7429x+19.8286，图示于图4中：x=12:16,1

20、8:2:26; y=52 55 58 61 65 70 75 80 85 91;% 定义透视电压和被透视得厚度 plot(x,y,'r+')%利用plot函数画出y随x变化的曲线hold on%保持原图形xlim(gca,6 28)%具有两个元素的数值向量xlabel('x'),ylabel('y')%xy轴的说明文字p= polyfit(x,y,1);% 得到最小二乘拟合多项式的系数 xx=10:0.2:28;%初值是10，步长是0.2，末值是28的行向量yy=polyval(p,xx)% 利用polyval函数计算函数近似值plot(xx,y

21、y) % 画出误差函数曲线   图四 y随x线性变化曲线 >> x=12:16,18:2:26; >> y=52 55 58 61 65 70 75 80 85 91; >> p=polyfit(x,y,1) %得到最小二乘拟合多项式的系数 p =      2.7429   19.8286 2、一元非线性经验公式  首先将初步选定的经验公式变换为直线式Y=AX+B式中，Y,X为只含有一个变量x或y的函数，A和B是与变换前经验公式

22、参数a、b有关的常数和系数。常用的变换如下面的表6： 表六：一元非线性经验公式序号非线性方程 线性化方程 线性化变量 1 2 3 下面我们通过一个具体的例子来说明求解的方法。例8.下面的测量数据（表7），若用指数形式的经验公式拟合，试计算公式的参数a与b。x 2 4 8 16 25 32 50 64 100 y 24.5 37 56.8 85.5 112.5 129.5 171.5 200 260.5 如表6中所示，用lny=lna+blnx对非线性方程线性化，写作 y=a+bx其中y=lny,x=lnx,a=lna即aMATLAB中程序如下生成图五的程序：x=2 4 8 16 25 32 50 64 100y=24.5 37 56.8 85.5 112.5 129.5 171.5 200 260.5;%给出两个矩阵xs=log(x);ys=log(y); % 将两个函数值赋值给两个变量p=polyfit(xs,ys,1); %得到最小二乘拟合多项式的系数a=exp(p(2);b=p(1)%定义函数的系数 plot(x,y,'r+')%用plot函数画图 hold on%保持原图像 xlabel('x'),ylabel('y')%x,y轴说明文字的句柄 xlim(gac,0