《数据分析与方程求解.ppt》由会员分享，可在线阅读，更多相关《数据分析与方程求解.ppt（34页珍藏版）》请在三一办公上搜索。
1、数据分析,主要内容,一、多项式计算二、数据的导入与导出三、数据的统计分析四、插值和拟合五、Mathematica插值六、Mathematica拟合,一、多项式计算,多项式可看做是符号表达式,利用符号运算可进行处理多项式的向量表示 设有多项式 表示为:(1)多项式的加法对应系数的加减运算(2)多项式的乘法：conv(P1,P2)(3)多项式除法：Q,r=deconv(P1,P2)其中Q返回相除的商，r返回余式。Q和r仍是多项式的系数向量(4)代数多项式求值：Y=polyval(P,x)若x为一数值，则求多项式P在该点的值；若x为向量或矩阵，则对向量或矩阵中的每个元素求多项式P的值(5)矩阵多项式
2、求值：Y=polyvalm(P,X)X为向方阵，结果为矩阵X的乘与和(6)多项式的求根：x=roots(P)在复数范围内求多项式P的全部根,二、数据的导入与导出,Matlab支持的数据文件Matlab特制的数据文件：.mat文件通用数据文件：文本文件、EXCEL文件、数据库文件 等特制数据文件的打开与保存保存右击需保存的变量，通过菜单操作save命令保存变量到数据文件打开通过菜单“File|open”load 命令通用数据文件通过“File|Import Data”打开导入向导,三、数据的统计分析,求最大值和最小值（max和min）Y,U=max(A)A为向量,将A的最大值存入U,最大值序号存
3、入UA为矩阵,Y表示每列的最大值,U记录每列最大值的行号 U=max(A,B)由A和B中对应元素的较大者组成矩阵U。求和、积、均值、中值、累加和、累乘sum和prod、mean和median、cumsum和cumprod 例：x=1,ones(1,10)*2;y=cumprod(x);sum(y)求得1+2+22+210排序Y,I=sort(A,dim,mode)Y是排序后的矩阵,I记录Y中的元素在A中的位置.dim取1(对列排序);取2(对行排序);默认列排序.mode的取值ascend(升序);descend(降序);默认升序.,向量和矩阵,四、插值与拟合,引例1:在一天24小时内,从零点
4、开始每间隔2小时测得的环境温度数据分别为(度)12,9,9,10,18,24,28,27,25,20,18,15,13.请推测中午1点（即13点）时的温度。引例2:在化学反应中,为研究某化合物的浓度随时间的变化规律,测得一组数据如下：请给出变化规律,由此推测t为1、1.5、2、2.5、10、10.5、11分钟时的值.利用已给出的值,计算相关的值设有一组实验观测数据(xi,yi),i=0,1,n,如何揭示自变量x与因变量y之间的关系？寻找近似的函数关系表达式y=f(x)常用的方法:插值与拟合,1、插值,设测得的n个点的数据为(xi,yi)(i=1,2,n),构造一个函数y=g(x),使得在xi(
5、i=1,2,n)有g(xi)=yi。g(x)称为插值函数。插值函数自变量的个数：一维插值、二维插值、多维插值等构造插值函数方法:线性插值、多项式插值、样条插值等一维数据的插值Y1=interpl(X,Y,X1,mothod)根据X和Y的值插值，并求插值函数在X1处的值Y是与X等长的向量,Y1是插值函数在X1处的值Y是矩阵,X与Y的每一列分别插值，并分别在X1的值X1是向量或标量,若X1中的元素不在X的范围内,则插值结果为NaNmethod的取值 nearest：最近插值、linear：线性插值 spline：三次样条插值、cubic：三次插值 二维数据插值 Z1=interp2(X,Y,Z,X
6、1,Y1,method)根据X、Y和Z的值插值，并求插值函数在X1、Y1处的值X、Y是向量(X:行,Y:列),Z为函数值(size(Z)=length(Y)length(X)X、Y、Z是矩阵,X与Y同维,Z(i,j)是在点(X(i,j)，Y(i,j)处的值,引例1的求解,在一天24小时内,从零点开始每间隔2小时测得的环境温度分别为(度)12,9,9,10,18,24,28,27,25,20,18,15,13.请推测中午1点（即13点）时的温度。解：(1)构造数据h=0:2:24;T=12,9,9,10,18,24,28,27,25,20,18,15,13;(2)了解数据的分布绘出数据的图形：p
7、lot(h,T,*)(3)选取适当插值方法，计算结果interp1(h,T,13,spline)interp1(h,T,13,cubic)插值函数的评价X=0:24;%选取更多的数据Y1=interp1(h,T,X,spline);Y2=interp1(h,T,X,cubic);subplot(2,1,1);plot(h,T,:,X,Y1);title(spline)subplot(2,1,2);plot(h,T,:,X,Y2,r);title(cubic),例：某实验对一根长为10米的钢轨进行热源的温度传播测试。如下表，其中x表示测量点，h表示测量时间，T表示测得的温度。试用线性插值求出在一
8、分钟内每隔20秒、钢轨每隔1米处的温度。解：x=0:2.5:10;h=0,30,60;T=95,14,0,0,0;88,48,32,12,6;67,64,54,48,41;x1=0:10;h1=0:20:60;T1=interp2(x,h,T,x1,h1)等价于X,H=meshgrid(0:2.5:10,0:30:60);T=95,14,0,0,0;88,48,32,12,6;67,64,54,48,41;X1,H1=meshgrid(0:10,0:20:60);T1=interp2(X,H,T,X1,H1),2、拟合,根据一组数据(xi,yi)(i=1,2,n)，要求确定一个函数y=f(x)
9、，使这些数据点与曲线y=f(x)总体来说尽量接近，称为曲线拟合。拟合的原理：最小二乘法 最小若拟合函数f(x)是一个多项式,称为多项式拟合.Matlab中多项式拟合命令：P=polyfit(X,Y,m)求数据X与Y的m阶拟合多项式P为拟合多项式的系数,长为m+1的向量X与Y是等长的向量通过polyval求拟合函数在自变量处的值,引例2的求解,引例2:在化学反应中,为研究某化合物的浓度随时间的变化规律,测得一组数据如下：请给出变化规律,推测t为1、1.5、2、2.5、10、10.5、11分钟时的值.解：(1)构造数据t=1:10;y=4,6.4,8,8.4,9.28,9.5,9.7,9.86,1
10、0,10.2;(2)选取拟合多项式的阶plot(t,y,o)(3)构造拟合多项式 p=polyfit(t,y,2)p=-0.1109 1.7922 2.9467(4)计算相关的值 t1=1:0.5:11;y1=polyval(p,t1);,拟合函数的评价选取不同的阶做拟合根据散点图直观判断计算最小二乘指标例t=1:10;y=4,6.4,8,8.4,9.28,9.5,9.7,9.86,10,10.2;p=polyfit(t,y,2);y2=polyval(p,t);plot(t,y,:o,t,y2,-*)sum(y-y2).2)类似做3次、4次多项式拟合，从中进行比较,例：在彩色显影中，由经验得
11、知，形成染料的光学密度与析出银的光学密度由公式确定，试验测得如下一批数据：求y关于x的拟合函数。解：由给定的经验公式来求拟合函数，无法直接用polyfit函数，我们换一种思路，将经验公式两边取对数，得 令Y=lny，X=1/x只需求Y关于X的线性拟合x=0.05,0.06,0.07,0.10,0.14,0.20,0.25,0.31,0.38,0.43,0.47;y=0.10,0.14,0.23,0.37,0.59,0.79,1.00,1.12,1.19,1.25,1.29;X=1./x;Y=log(y)P=polyfit(X,Y,1)P=-0.1459 0.5476则Y=-0.1459X+0.
12、5476y=exp(Y)=exp(0.5476)e(-0.1459X)=exp(0.5476)e(-0.1459/x),整理即可,五、Mathematica插值,构造插值对象的函数：Interpolationdata,InterpolationOrdern功能：对数据data进行插值,并可设置插值多项式的次数n,默认值为3。注:生成一个InterpolatingFunction插值范围，目标所得目标为近似函数,不显示所构造的函数。数据表示：x0,f0,x1,f1,xn,fn(平面点的坐标)引例1求解解：data=0,12,2,9,4,9,6,10,8,18,10,24,12,28,14,27,
13、16,25,18,20,20,18,22,15,24,13ListPlotdataf=Interpolationdataf(13),六、Mathematica拟合,格式：Fit数据，拟合函数的基，变量功能：用数据data，按给定的变量和拟合函数的基构造拟合函数。常用的几种格式：Fitdata,1,x,x用数据data作线性拟合函数a+bxFitdata,Tablexi,i,0,n,x作n次多项式拟合引例2的求解：构建数据表：data=1,4,2,6.4,3,8.0,4,8.4,5,9.28,6,9.5,7,9.7,8,9.86,9,10,10,10.2选取拟合的基。画图，分别作3次、4次、5次
14、、6次多项式拟合，求出其误差，比较、判断几次拟合较为合适。,插值与拟合插值要求函数在每一个观测点处一定要满足yi=f(xi).拟合主要考虑到观测数据受随机误差的影响,寻求整体误差最小、较好反映观测数据的近似函数.插值与拟合一般过程：构造数据确定插值多项式的阶（可采用经验公式或画图）选取命令插值求出并分析所得结果数据拟合的一般过程：考察数据来源,确定有无经验公式；若有经验公式,则可直接进行拟合（必要时可对数据作些变换处理）;若无经验公式,选择拟合函数的类型（画图）求出拟合函数从数据、图形等途径对拟合结果分析;确定结果或调整拟合函数,方程求解,主要内容,符号方程的求解 符号代数方程求解 符号常微分
15、方程求解 线性方程组求解 利用左除运算符 利用矩阵的分解 非线性方程求解 单个自变量的非线性方程求解 非线性方程组求解 微分方程求解 Mathematica方程的求解,一、符号方程的求解,符号代数方程求解 代数方程是指未涉及微积分运算的方程 求解命令solve(s)求解符号代数方程，求解变量是默认变量。solve(s1,s2,.,sn,v1,v2,.,vn)求解符号表达式s1,s2,.,sn所组成的代数方程组的解，求解变量分别为v1,v2,.,vn。例：solve(x2-a*x-b=0)solve(x2-a*x-b)solve(x2-a*x-b=0,x)syms x a b f=x2-a*x-
16、b;solve(f,x);solve(f)syms x y z x y z=solve(x+2*y-z=27,x+z=3,x2+3*y2=12,x,y,z)f=x5-a*x-b=0;solve(f,x)%无法求出解析解 f=x5-x-3;solve(f)%求近似解,符号常微分方程求解方程的表示：Dy表示y，D2y表示y 求解命令dsolve(e,c,v)求解微分方程e在初始条件c下的特解,参数v描述方程中的自变量.例：(1)求 当y(0)=1时的特解。y=dsolve(Dy=2*y*y2,y(0)=1,x)(2)求 的通解。x,y=dsolve(Dx=4*x-2*y,Dy=2*x-y,t),二
17、、线性方程组求解,若有线性方程组利用左除运算符A=2,1,-1,1;1,-5,0,7;0,2,1,-1;1,6,-1,-4;b=13,-9,6,0;x=Ab 利用矩阵的分解LU分解：将矩阵表示为一个下三角矩阵与一个上三角矩阵的乘积。L,U=lu(X):产生L和U,使得X=LU。A=2,1,-1,1;1,-5,0,7;0,2,1,-1;1,6,-1,-4;b=13,-9,6,0;L,U=lu(A);x=U(Lb)QR分解:是将矩阵分解为一个正交矩阵和一个上三角矩阵的乘积.：Q,R=qr(X):产生Q和R,使得X=QR。,三、非线性方程求解(近似解),例：考察函数作图观察其根的分布尝试用solve
18、求解单个自变量的非线性方程求解 z=fzero(fun,x0)z=fzero(fname,x0)在x0附近寻找函数fun的近似根(fname是待求根的函数名)求上述函数在区间-3,1内的根（1）画图确定根的范围 f=abs(x*sin(x)-exp(x)-1;fplot(f,-3,1)观察图形得出在x=0、-1.2、-2.6附近有根（2）调用fzero函数分别求出其根 fzero(f,0)fzero(f,-1.2)fzero(f,-2.6),上述问题也可这样求解(1)定义函数文件 funx.m,function y=funx(x)y=abs(x.*sin(x)-exp(x)-1;(2)画图确定
19、根的位置 fplot(funx,-3,1)%funx为函数funx的句炳观察图形得到出在x=0、-1.2、-2.6附近有根(3)调用fzero函数分别求出其根 fzero(funx,0)fzero(funx,-1.2)fzero(funx,-2.6)上述问题还可这样求解(1)定义内联（匿名）函数ff=(x)abs(x.*sin(x)-exp(x)-1;(2)画图确定根的位置 fplot(ff,-3,1)(3)调用fzero函数分别求出其根 fzero(ff,-2.5)fzero(ff,-1.5)fzero(ff,0),非线性方程组求解,求解方程组 初始值为-5,5MATLAB的优化工具箱（Op
20、timization Toolbox）对于非线性方程组F(x)=0,求解命令为：X=fsolve(fname,X0,options)解(1)定义函数文件myfun.m function f=myFun(x)f(1)=2*x(1)-x(2)-exp(-x(1);f(2)=-x(1)+2*x(2)-exp(-x(2);(2)提交给求解方程的函数-5,-5%初始值 op=optimset(display,off);%修改参数选项,中间结果不显示 x=fsolve(myfun,x0,op)x=0.5671 0.5671,四、微分方程求解,求解微分方程：MATLAB提供了多个求常微分方程数值解的函数,格
21、式为:t,y=solver(fname,tspan,y0)solver为求常微分方程数值解的函数t和y分别给出时间向量和相应的状态向量。fname是定义f(t,y)的函数文件名，该函数文件必须返回一个列向量。tspan为求解区间形式为t0,tf,y0是初始状态列向量。解：建立函数文件funt.mfunction yp=funt(t,y)yp=(y2-t-2)/4/(t+1);求解微分方程t,y=ode23(funt,0,10,2);%求数值解,五、Mathematica方程的求解,方程的表示方程中的等号（=）应用逻辑等号（=）根的表示：逻辑量、转换规则多项式方程求解Solve方程或方程组,变量
22、列表NSolve方程或方程组,变量列表注：Solve主要是处理多项式方程，尽可能给出精确解，若给不出显示公式解，给出隐式公式解。例：求方程x5-x+11=0的根。NSolve求出的是近似数值解(可求任意多项式方程)方程组的表示：方式1：lhs1=rhs1,lhs2=rhs2方式2：lhs1=rhs1&lhs2=rhs2&方式3：lhs1,lhs2,=rhs1,rhs2.,Mathematica方程的求解,一般方程的求解FindRoot方程,变量,初值FindRoot方程,变量,min,max在初值附近或区间内找方程的根，一次只能算出一个近似数值解。求y=|x*sin(x)-ex|-1在区间-3
23、,1内的根。,最优化问题求解,主要内容,一、无约束问题求解（求极值、求最值）模型:minxf(x)求一组x(x=x1,x2,xnT)使得目标函数f(x)为最小二、有约束问题求解模型：求一组x(x=x1,x2,xnT)使得目标函数f(x)为最小,且x满足约束条件G(x)0.约束条件可表示为：线性不等式约束：Ax b线性等式约束：Aeqx=beq非线性不等式约束：Cx b非线性等式约束：Ceqx=0X的上界和下界：Lbnd x Ubng 三、线性规划问题求解线性约束条件下线性目标函数的极值问题,一、无约束问题求解（求极值、最值）,若函数可导,换转为求导数为零的点。例:求y=x*sin(x)-ex在
24、区间-3,2内的极值点。解:fplot(x*sin(x)-exp(x),-3,2)%图形显示在x=-2附近有一极大值点 syms x f=x*sin(x)-exp(x);g=diff(f)%求导函数 sin(x)+x*cos(x)-exp(x)fzero(sin(x)+x*cos(x)-exp(x),-2)结果：-2.0745故在x-2.0745处函数取得极大值。上述问题也求解如下：syms x f=x*sin(x)-exp(x);fplot(char(f),-3,2);%将符号表达式转换为字符串 x0=fzero(char(diff(f),-2);y0=subs(f,x,x0);%将符号表达
25、式f中的x替换为x0,求得极大值.x0,y0-2.0745,1.6912,Matlab专门提供了求极小值的函数，格式为：(1)x,fval=fminbnd(fname,x1,x2,option)求一元函数在区间(x1,x2)中的极小值点x和极小值。(2)x,fval=fminsearch(fname,x0,option)用单纯形法求多元函数在x0附近的极小值点x和极小值.(3)x,fval=fminunc(fname,x0,option)用拟牛顿法求多元函数在x0附近的极小值点x和极小值.例:求y=x*sin(x)-ex在区间-3,2内的极值点。(1)建立函数文件myfin.m,命令如下：fu
26、nction y=myfin(x)y=-(x.*sin(x)-exp(x);(2)调用fminbnd函数求极大值,命令如下：x,fval=fminbnd(myfin,-3,-1.5)结果：x=-2.0745，fval=-1.6912y=x*sin(x)-ex在x=-2.0745处取得极大值1.6912。等价于f=(x)-(x*sin(x)-exp(x);x,fval=fminbnd(f,-3,-1.5),二、有约束问题求解,模型：求一组x(x=x1,x2,xnT)使得目标函数f(x)为最小,且x满足约束条件G(x)0.约束条件可表示为：线性不等式约束：Ax b线性等式约束：Aeqx=beq非线
27、性不等式约束：Cx b,非线性等式约束：Ceqx=0X的上界和下界：Lbnd x Ubng Matlab优化工具箱提供了解各种约束下的最优化问题。格式：x,fval=fmincon(fname,x0,A,b,Aeq,beq,Lbnd,Ubnd,NonF,options)求目标函数在x0附近满足约束条件的极小值点和极小值A,b：构成线性不等式约束Aeq,beq：线性等式约束Lbnd,Ubnd：x的上下界常用：x,fval=fmincon(fname,x0,A,b)x,fval=fmincon(fname,x0,A,b,Lbnd,Ubnd),例:设求解有约束最优化问题：解:(1)编写目标函数文件f
28、op.mat,命令如下：function f=fop(x)f=0.4*x(2)+x(1)2+x(2)2-x(1)*x(2)+1/30*x(1)3;(2)设定约束条件，并调用fmincon求解。x0=0.5;0.5;%列向量A=-1,-0.5;-0.5,-1;b=-0.4;-0.5;lb=0;0;%列向量options=optimset(Display,off);x,f=fmincon(fop,x0,A,b,lb,options),三、线性规划问题求解,线性规划问题是研究线性约束条件下线性目标函数极值问题的数学理论和方法。线性规划问题的标准形式为：在Matlab中求解线性规划问题使用函数linprog，格式为 x,fval=linprog(f,A,b,Aeq,beq,lb,ub)x为最优解，fval为目标函数的最优值，f为目标函数的系数向量。例：求解线性规划问题：解：f=2;1;A=-3,-1;-4,-3;-1,-2;b=-3;-6;-2;lb=0;0;x,fval=linprog(f,A,b,lb)x=0.6000 1.2000fval=2.4000,
