第四讲数据分析方法

1、第四讲 数据分析方法第一节、数据拟合问题：给定一批数据点（输入变量与输出变量的数据），需确定满足特定要求的曲线或曲面。如果输入变量和输出变量都只有一个，则属于一元函数的拟合和插值；而若输入变量有多个，则为多元函数的拟合和插值（有点回归分析的意思）解决方案：（1） 若要求所求曲线（面）通过所给所有数据点，就是插值问题；（2） 若不要求曲线（面）通过所有数据点，而是要求它反映对象整体的变化趋势，这就是数据拟合，又称曲线拟合或曲面拟合。注意：插值和拟合都是要根据一组数据构造一个函数作为近似，由于近似的要求不同，二者的数学方法上是完全不同的。而面对一个实际问题，究竟应该用插值还是拟合，有时容易确定，有

2、时则并不明显。例 1：下面数据是某次实验所得，希望得到 X和 f之间的关系？xf12479121315171.53.96.611.715.618.819.620.621.1曲线拟合问题最常用的解法最小二乘法的基本思路第一步：确定拟合的函数类型 y = f (x;a1,a2,am )，其中a1,a2,am为待定系数。（函数类型的确定可以根据内在的规律确定，如果无现成的规则，则可以通过散点图，联系曲线的形状进行分析）第二步：确定a1,a2,am的最小二乘准则：要求n个已知点(xi, yi )与曲线 y = f (x)n的距离di的平方和 (yi f x 2 最小 。( )ii=1用 MATLAB作

3、拟合1多项式拟合。作多项式 y = a0xm+a1x m1 + am拟合，可利用a=polyfit(x,y,m)其中 x,y为给出的数据，m为多项式的次数。多项式在 x处的值 y可用以下命令计算：y=polyval（a,x）2用 MATLAB作非线性最小二乘拟合Matlab的提供了两个求非线性最小二乘拟合的函数：lsqcurvefit和 lsqnonlin。两个命令都要先建立 M-文件 fun.m，在其中定义函数 f(x)。（1）x = lsqcurvefit (fun,x0,xdata,ydata);（2）x =lsqcurvefit (fun,x0,xdata,ydata,options)

4、;（3）x = lsqcurvefit (fun,x0,xdata,ydata,options,grad);（4）x, options = lsqcurvefit (fun,x0,xdata,ydata,);（5）x, options,funval = lsqcurvefit (fun,x0,xdata,ydata,);1（6）x, options,funval, Jacob = lsqcurvefit (fun,x0,xdata,ydata,);（7）x=lsqnonlin（fun，x0）；（8）x= lsqnonlin （fun，x0，options）；（9）x= lsqnonlin （f

5、un，x0，options，grad）；（10）x，options= lsqnonlin （fun，x0，）；（11）x，options，funval= lsqnonlin （fun，x0，）；在使用 lsqcurvefit与 lsqnonlin命令时，共同的问题是要先知道函数的类型，而拟合其实是决定函数中的待定系数。第二节 插值插值的基本问题：给出n个数对，(Pi, f (Pi),i =1,2,n，求点 P处对应的函数值 f (P)。一、一维插值已知 n +1个节点 (xi, yi ),i = 0,1,n，求任意点 *处的函数值 *。常用的插值方法xy有拉格朗日多项式插值、牛顿插值、分段线性

6、插值、Hermite插值和三次样条插值。9 分段线性插值：将各数据点用折线连接起来9 多项式插值：求一个多项式通过所有数据点，可以假设出多项式的系数，最后通过求解方程得到每个系数（拉格朗日插值，用n次多项式描述 n +1个点）9 样条插值：分段多项式的光滑连接（三次样条插值）9 牛顿插值：利用节点之间的各阶差商和差分构造多项式9Hermite插值：对插值函数，不仅要求它在节点处与函数同值，而且要求它与函数有相同的一阶、二阶甚至更高阶的导数值（1）MATLAB命令：y=interp1(x0,y0,x,method)method指定插值的方法，默认为线性插值。其值可为：nearest 最近项插值l

7、inearsplinecubic线性插值立方样条插值立方插值。所有的插值方法要求 x0是单调的。当 x0为等距时可以用快速插值法，使用快速插值法的格式为*nearest、*linear、*spline、*cubic。（2）三次样条插值在 Matlab中的实现在 Matlab中数据点称之为断点。如果三次样条插值没有边界条件，最常用的方法，就是采用非扭结（not-a-knot）条件。这个条件强迫第 1个和第 2个三次多项式的三阶导数相等。对最后一个和倒数第 2个三次多项式也做同样地处理。Matlab中三次样条插值也有现成的函数：y=interp1(x0,y0,x,spline)；y=spline(

8、x0,y0,x)；pp=csape(x0,y0,conds),pp=csape(x0,y0,conds,valconds)，y=ppval(pp,x)。其中 x0,y0是已知数据点，x是插值点，y是插值点的函数值。对于三次样条插值，我们提倡使用函数 csape，csape的返回值是 pp形式，要求插值点2的函数值，必须调用函数 ppval。pp=csape(x0,y0)：使用默认的边界条件，即 Lagrange边界条件。pp=csape(x0,y0,conds,valconds)中的 conds指定插值的边界条件，其值可为：complete边界为一阶导数，一阶导数的值在 valconds参数中

9、给出，若忽略 valconds参数，则按缺省情况处理。not-a-knot 非扭结条件periodicsecond周期条件边界为二阶导数，二阶导数的值在 valconds参数中给出，若忽略 valconds参数，二阶导数的缺省值为0, 0。variational 设置边界的二阶导数值为0,0。对于一些特殊的边界条件，可以通过 conds的一个1 2矩阵来表示，conds元素的取值为 0，1，2。conds(i)=j的含义是给定端点i的 j阶导数，即 conds的第一个元素表示左边界的条件，第二个元素表示右边界的条件，conds=2,1表示左边界是二阶导数，右边界是一阶导数，对应的值由 valc

10、onds给出。例 2 机床加工待加工零件的外形根据工艺要求由一组数据(x, y)给出（在平面情况下），用程控铣床加工时每一刀只能沿 x方向和 y方向走非常小的一步，这就需要从已知数据得到加工所要求的步长很小的(x, y)坐标。表中给出的 x, y数据位于机翼断面的下轮廓线上，假设需要得到 x坐标每改变 0.1时的y坐标。试完成加工所需数据，画出曲线，并求出 x = 0处的曲线斜率和13 x 15范围内 y的最小值。xy03579112.0121.8131.2141.0151.601.21.72.02.1要求用分段线性和三次样条计算。解 编写以下程序：x0=0 3 5 7 9 11 12 13

11、14 15;y0=0 1.2 1.7 2.0 2.1 2.0 1.8 1.2 1.0 1.6;x=0:0.1:15;y2=interp1(x0,y0,x);y3=interp1(x0,y0,x,spline);pp1=csape(x0,y0);y4=ppval(pp1,x);pp2=csape(x0,y0,second);y5=ppval(pp2,x);x,y1,y2,y3,y4,y5subplot(2,2,1)plot(x0,y0,+,x,y2)title(Piecewise linear)subplot(2,2,2)plot(x0,y0,+,x,y3)3title(Spline1)subp

12、lot(2,2,3)plot(x0,y0,+,x,y4)title(Spline2)dx=diff(x);dy=diff(y3);dy_dx=dy./dx;dy_dx0=dy_dx(1)ytemp=y3(131:151);ymin=min(ytemp);index=find(y3=ymin);xmin=x(index);xmin,ymin二、二维插值：前面讲述的都是一维插值，即节点为一维变量，插值函数是一元函数（曲线）。若节点是二维的，插值函数就是二元函数，即曲面。如在某区域测量了若干点（节点）的高程（节点值），为了画出较精确的等高线图，就要先插入更多的点（插值点），计算这些点的高程（插值）。

13、（1）（2）网格节点：已知m n个节点 (xi , y j , zij ),i = 1,2,m; j = 1,2,n，构造一个二元函数 z = f (x, y)通过全部已知节点，即 zij = f (xi , y j )，再利用z = f (x, y)插值。散乱节点：已知 n 个节点 (xi, yi, zi ),i = 1,2,n ，构造一个二元函数z = f (x, y)通过全部已知节点，即 zi = f (xi, yi )，再利用 z = f (x, y)插值。二维插值方法：（1）最邻近插值：二维或高维情形的最邻近插值，与被插值点最邻近的节点的函数值即为所求。（2）分片线性插值：将四个插值

14、点（矩形的四个顶点）处的函数值依次简记为：f (xi, yj)=f1，f (xi+1, yj)=f2，f (xi+1, yj+1)=f3，f (xi, yj+1)=f4分两片的函数表达式如下：第一片(下三角形区域)： f (x, y) = f1 + ( f 2 f1)(x xi ) + ( f3 f 2)(y y j )第二片(上三角形区域)： f (x, y) = f1 + ( f 4 f1)(y y j ) + ( f3 f 4)(x xi )（3）双线性插值：双线性插值是一片一片的空间二次曲面构成。双线性插值函数的形式如下： f (x, y) = (ax + b)(cy + d)。其中有

15、四个待定系数，利用该函数在矩形的四个顶点（插值节点）的函数值，得到四个代数方程，正好确定四个系数。用 MATLAB作网格节点数据的插值：4z=interp2(x0,y0,z0,x,y,method) method 可以选择 nearest：最邻近插值；linear：双线性插值；cubic：双三次插值；缺省时, 双线性插值。如果是三次样条插值，可以使用命令pp=csape(x0,y0,z0,conds,valconds)，z=fnval(pp,x,y)其中 x0,y0分别为m维和n维向量，z0为m n维矩阵， z为矩阵，它的行数为 x的维数，列数为 y的维数，表示得到的插值，具体使用方法同一维插

16、值。例3 在一丘陵地带测量高程， x和 y方向每隔100米测一个点，得高程如下表，试插值一曲面，确定合适的模型，并由此找出最高点和该点的高程。xy100200300400100200697712674626300624630598552400 500478 450478 420412 400334 310636698680662解 编写程序如下：x=100:100:500;y=100:100:400;z=636 697 624 478 450;698 712 630 478 420;680 674 598 412 400;662 626 552 334 310;pp=csape(x,y,z)x

17、i=100:10:500;yi=100:10:400cz1=fnval(pp,xi,yi)cz2=interp2(x,y,z,xi,yi,spline)i,j=find(cz1=max(max(cz1)x=xi(i),y=yi(j),zmax=cz1(i,j)用 MATLAB作散点数据的插值计算cz =griddata（x,y,z,cx,cy,method）method可以选择nearest：最邻近插值；linear ：双线性插值；cubic：双三次插值；v4：Matlab提供的插值方法；缺省时, 双线性插值例 4 在某海域测得一些点（x,y）处的水深 z由下表给出，在矩形区域（75,200）

18、(-50,150)内画出海底曲面的图形。x 129y 7.5z 4140141.58103.5236881478185.522.56195137.58105157.56.59107.5-819778156.5816266.59162117.5-33.5985.5384884解 编写程序如下：x=129,140,103.5,88,185.5,195,105,157.5,107.5,77,81,162,162,117.5;y=7.5,141.5,23,147,22.5,137.5,85.5,-6.5,-81,3,56.5,-66.5,84,-33.5;z=4,8,6,8,6,8,8,9,9,8,8

19、,9,4,9;xi=75:1:200;5yi=-50:1:150;zi=griddata(x,y,z,xi,yi,cubic)subplot(1,2,1)plot(x,y,*)subplot(1,2,2)mesh(xi,yi,zi)第三节、回归分析方法“回归”问题最早来源于生物界，英国生物学家兼统计学家高尔顿（Galton,1822-1911）发现同一种族中儿子的平均高度介于其父亲的高度与种族平均高度之间。儿子的身高有返归于种族平均身高的趋势，即回归于种族的平均身高。回归分析是指对具有相关关系的现象，根据其关系形态，选择一个合适的数学模型，用来近似地表示变量间的平均变化关系的一种统计方法。回归

20、分析的分类：按照回归模型中变量个数分（一元回归，多元回归）；按照回归曲线的形态分（线性回归，非线性回归）；按照是否要求总体分布类型已知分（参数回归，非参数回归）一元线性回归特点 两个变量中，一个是自变量，一个是因变量 回归方程不是抽象的数学模型，而随机方程，可以进行实证 因果关系不明显时，应同时作两个回归方程 回归系数具有较强的经济含义 作为回归模型的因变量是随机变量，而自变量是确定性变量，即可控变量1一元线性回归模型一般地，称由 y = 0 + 1x + 确定的模型为一元线性回归模型，记为y = 0 + 1x +E = 0, D = 2其中固定的未知参数 0， 1称为回归系数，自变量 x也称

21、为回归变量， Y = 0 + 1x称为 y 对 x的回归直线方程。一元线性回归分析的主要任务1用试验值（样本值）对 0， 1和 作点估计；2对回归系数 0， 1作假设检验；3在 x = x0处对 y作预测，对 y作区间估计。回归系数的最小二乘估计有n组独立观测值，(xi, yi ),i = 1,2,n。yi = 0 + x1 + i,i = 1,2,.,n设 E i = 0, D i = 2 且1 2，., n相互独立6nn = (yi 0 1xi)2 。最小二乘法就是选择2i记Q = Q( 0,1) =0和 1的估i=1i=1计 0， 1使得Q(0,1) = min Q( 0,1)。计算得到

22、 ,01 0 = y 1xxy xy =1x x2 21n1n1n1n其中 x =xi, y = nyi ， x 2 =xi 2, xy =xi yi 。（经验）回归方程为：nnni=1i=1i=1i=1y = 0 + 1x = y + 1(x x) 。记 Q = Q(0, )= (yi 0 1xi)2 = (yi y 2 ，称 e为残差平方和或剩nnQ )ie1i=1i=1余平方和。 2的无偏估计为 e2= Qe (n 2) 。回归方程的显著性检验对回归方程 Y = 0 + 1x的显著性检验，归结为对假设 H 0 : 1 = 0;H1 : 1 0进行检验。假设 H 0 : 1 = 0被拒绝，

23、则回归显著，认为 y与 x存在线性关系，所求的线性回归方程有意义；否则回归不显著，y与 x的关系不能用一元线性回归模型来描述，所得的回归方程也无意义。UQe /(n 2)nF检验法：当 H 0成立时， F = F(1,n 2)，其中 U = (yi y)2（回i=1归平方和）。若 F F1 (1,n 2)，拒绝 H 0，否则就接受 H 0。回归系数的置信区间0和 1置信水平为1 的置信区间分别为1+ x ,0 + t (n 2) e21+ x 2 0 t (n 2) e和nLxxnLxx1122 t (n 2) e / Lxx ,1 + t (n 2) e / L ；11212xx7QeQe2

24、的置信水平为1 的置信区间为 2 (n 2) , (n 2)。21 22预测用 y 的回归值 y0 = 0 + 1x00作为 0的预测值。 0的置信水平为yy1的预测区间为(x0 x)Lxx2y0 (x0), y0 + (x0)，其中 (x0) = et (n 2) 1+ 1 +。特别，当 n很1n2大 且 x0 在 x 附 近 取 值 时 ,y 的 置 信 水 平 为 1 的 预 测 区 间 近 似 为y eu , y + eu。12122可线性化的一元非线性回归（需要配曲线）先对两个变量 x和 y作n次试验观察得(xi , yi ),i =1,2,.,n画出散点图，根据散点图确定须配曲线的

25、类型.然后由n对试验数据确定每一类曲线的未知参数a和b。采用的方法是通过变量代换把非线性回归化成线性回归，即采用非线性回归线性化的方法。通常选择的六类曲线如下（1）双曲线 1 = a +byx（2）幂函数曲线 y = axb，其中 x 0,a 0（3）指数曲线 y = aebx其中参数 a 0（4）倒指数曲线 y = aeb / x其中 a 0（5）对数曲线 y = a + bln x, x 01（6）S型曲线 y = a + bex3多元线性回归Y = X +一般称 为高斯马尔柯夫线性模型(k 元线性回归模E() = 0,COV(,) = 2I n型)，并简记为(Y, X, 2I n )。8

26、 y1 1 x11 x12 . x1k 0 1 .1x21 x22 . x2k其 中 Y = ，X = ， =， = ，12 . . . . . . . . y1xn1 xn2 . xnkn nky = 0 + 1x1 +.+ kxk称为回归平面方程。线性模型(Y, X,2I n )考虑的主要问题是：（1）用试验值（样本值）对未知参数 和作点估计和假设检验，从而建立 y与 x1, x2,., xk2之间的数量关系；（2）在 x1 = x01, x2 = x02,., xk = x0k处对 y的值作预测与控制，即对 y作区间估计。对 i和2作估计，用最小二乘法求 0,., k 的估计量：作离差平

27、方和nQ =(yi 0 1xi1 . k xik )2i=1选择 0,., k 使Q达到最小。解得估计值 = (X T X)1(X TY)，得到的 代入回归平面i方程得 y = 0 + 1x1 + .+ k xk ，称为经验回归平面方程. i称为经验回归系数4多元非线性回归多项式回归设变量 x,Y 的回归模型为Y = 0 + 1x + 2x其中 p 是已知的， i (i = 1,2, p) 是未知参数， 服从正态分布 N(0, 2) 。Y = 0 + 1x + 2x + .+ x+ .+ p x + 2 p2k称为回归多项式，上面的回归模型称为多项式回归。k令 xi = xi，i=1，2，k多

28、项式回归模型变为多元线性回归模型.多元线性回归中的检验与预测假设 H 0 : 0 = 1 = . = k = 0。U / k当 H 0成立时， F = F(k,n k 1)。如果 F F1 (k,n k 1)，Qe /(n k 1)则拒绝 H 0，认为 y 与 x1, xk之间显著地有线性关系；否则就接受 H 0，认为 y 与 x1, xk之间线性关系不显著。9求 出 回 归 方 程 y = 0 + 1x1 +.+ k xk ， 对 于 给 定 自 变 量 的 值x*1,., xk， 用= 0 + 1x1+.+ k xk来预测 y= 0 + 1x1+.+ k xk+y.称 *为 *的点预测。yy*5逐步回归分析“最优”的回归方程就是包含所有对Y 有影响的变量, 而不包含对Y 影响不显著的变量回归方程。选择“最优”的回归方程有以下几种方法：（1）从所有可能的因子（变量）组合的回归方程中选择最优者；（2）从包含全部变量的回归方程中逐次剔除不显著因子；（3）从一个变量开始，把变量逐个引入方程；（4）“有进有出”的逐步回归分析逐步回归分析法的思想：（1）从一个自变量开始，视