数学建模第四章-1

4.1 预备知识：极值
4.1 预备知识：最小二乘拟合

数据：男子110米栏记录
成绩 13秒24 13秒21 13秒16 13秒00 12秒93 12秒92 12秒91 选手 米尔布恩 卡萨那斯 内赫米亚 内赫米亚 内赫米亚 金多姆 杰克逊 国籍 美国 古巴 美国 美国 美国 美国 英国 日期 1972年9月7日 1977年8月21日 1979年4月14日 1979年5月6日 1981年8月19日 1989年8月16日 1993年8月20日
>> [x,f,h]=fsolve(fun,-1.6) >> [x,f,h]=fsolve(fun,-0.6)
例4 求方程组在原点附近的解

使用 函数 句柄
1 x 4x y e 1 x x(1) 10 y x(2) 1 2 x 4 y x 0 8 function f=eg4_4fun(x) f(1)=4*x(1)-x(2)+exp(x(1))/10-1; f(2)=-x(1)+4*x(2)+x(1)^2/8;
牛顿法程序newton.m
function x=newton(fname,dfname,x0,e) if nargin<4, e=1e-4;end %默认精度 x=x0;x0=x+2*e; %使while成立，且进入 while后x0得到赋值 while abs(x0-x)>e x0=x;
x=x0-fname(x0)/dfname(x0);
[x,f]=fminsearch(fun,x0) x返回多元函数fun在初始值x0 附近的局部极小值点, f返回局部极小值. x, x0均为向量。
例 5 .求二元函数f(x,y)= 5-x4-y4+4xy在 原点附近的极大值。 解：max fmin(-f) x x(1), y x(2)
>> fun=inline('x(1)^4+x(2)^4-4*x(1)*x(2)-5'); >> [x,g]=fminsearch(fun,[0,0])
4.2 函数零点MATLAB指令
x=fzero(Fun, x0) 返回一元函数Fun的一个零点. x0为标量时, x 返回函数在x0附近的零点； x0为向量[a, b]时, x 返回在[a,b]中的零点(两端异号)
[x,f,h]=fsolve(Fun, x0) x: 返回多元函数Fun在x0附近的一个零 点，其中x, x0均为向量; f: 返回Fun在零点的函数值, 应该接近0; h: 返回值如果大于零，说明计算结果可 靠，否则计算结果不可靠。

用lsqcurvefit解例4.2
>>fun2 = inline('c(1)*x.^2+c(2)*x+c(3)','c','x') >> x=[0.1, 0.2, 0.15, 0, -0.2, 0.3]; >> y=[0.95,0.84,0.86,1.06,1.50,0.72];
>> c=lsqcurvefit(fun2,[0,0,0],x,y)
4.1 预备知识：最小二乘拟合
• 假设已知经验公式y=f(c,x)(c为参数, x为自变量), 要求 根据一批有误差的数据(xi,yi), i=0,1,…,n, 确定参数c. 这样的问题称为数据拟合。 • 最小二乘法就是求c使得均方误差最小化 Q(c)=
2 ( y f ( c , x )) i i i 0 n
end
例6 求方程 x 2 - 3 x + e x = 2 的正根 (要求精度 = 10 -6) 解 令f (x) = x 2 - 3 x + e x - 2, f(0)=-1, 当x > 2, f (x) > 0, f ’(x) > 0 即f (x)单调上升，所以根在[0,2]内。 先用图解法找初值， 再用牛顿法程序newton.m求解。
f (x) = 0， x=(x1, x2, …, xn), f=(f1, f2, …, fm)
4.1 预备知识：极值
设x为标量或向量，y=f(x)是xD上的标量值函数。 如果对于包含x=a的某个邻域 ,有 f(a)f(x) (f(a)f(x))对任意x成立, 则称a为f(x)的一个局 部极小(大)值点。 如果对任意xD，有f(a)f(x)（f(a)f(x)）成立， 则称a为f(x)在区域D上的一个全局极小(大)值点。
>> [x,f,h]=fsolve(@eg4_4fun,[0 0])
%使用Inline函数
fun=inline('[4*x(1)-x(2)+ exp(x(1)) /10 - 1, -x(1) + 4*x(2)+x(1)^2/8]','x');
[x,f,h]=fsolve(fun,[0 0])
%使用匿名函数
MATLAB中一个多项式用系数降幂排 列向量来表示。
例1.求多项式x3 + 2 x2 - 5的根 »p=[1 2 0 -5]; x=roots(p) , polyval(p,x) 例2.用2次多项式拟合下列数据. x 0.1 0.2 0.15 0 -0.2 0.3 y 0.95 0.84 0.86 1.06 1.50 0.72 M文件eg4_2.m
• 当f关于c是线性函数,问题转化为一个线性方程组求解。 • 如果f关于c是非线性函数，问题转化为函数极值问题
4.2 函数零点MATLAB指令
• 多项式
y=polyval(p,x) 求得多项式p在x处的值y，x可以是一个或多 个点(polynomial value) p3=conv(p1,p2) 返回多项式p1和p2的乘积(Convolution and polynomial multiplication) [p3,r]=deconv(p1,p2) p3返回多项式p1除以p2的商，r返回余 项 x=roots(p) 求得多项式p的所有复根. p=polyfit(x,y,k)用k次多项式拟合向量数据(x, y),返回多项式 的降幂系数
注：在使用fsolve, fminsearch等指令时， 多变量必须合写成一个向量变量，如用 x(1), x(2),…。
4.2 最小二乘拟合MATLAB指令
• 假设已知经验公式y=f(c,x)(c为参数, x为自变量), 要求 根据一批有误差的数据(xi,yi), i=0,1,…,n, 确定参数c. 这样的问题称为数据拟合。 • 最小二乘法就是求c使得均方误差最小化 Q(c)=
12秒88 12秒87
刘翔 罗伯斯
中国 古巴
2006年7月12日 2008年6月12日

问题：估计什么时候突破12秒50？
4.1 预备知识：最小二乘拟合
13.3 13.25 13.2 13.15 13.1 13.05 13 12.95 12.9 12.85 1970 1975 1980 1985 1990 1995 2000 2005 2010
fun=@(x)[4*x(1)-x(2)+exp(x(1))/10-1, -x(1) +4*x(2)+x(1)^2/8]; [x,f,h]=fsolve(fun,[0 0])
4.2 函数极值MATLAB指令
min(y) max(y) 返回向量y的最小值 返回向量y的最大值
[x,f]=fminbnd(fun,a,b) x返回一元函数fun在[a,b]内的 局部极小值点, f返回局部极小值
4.2 非线性函数的MATLAB表达
Fun=@Mfun 定义一个函数句柄，这里 Mfun是函数的M文件表达方式
Fun=inline(‘funstr’,’var’) 定义一个Inline 函数,其中funstr是函数的表达式, var是变量名
Fun=@(var)funstr 定义匿名函数，其中var 是变量名, funstr是函数的表达式
c=lsqcurvefit(Fun2,c0, x, y) 从外部输入数据， 这里Fun2为两变量c和x的函数 f(c, x)
用lsqnonlin解例4.2
先写M函数fitf.m function e=fitf(c) x=[0.1 0.2 0.15 0 -0.2 0.3]; y=[0.95 0.84 0.86 1.06 1.50 0.72]; e=y-(c(1)*x.^2+c(2)*x+c(3)); 然后执行 >> c=lsqnonlin(@fitf,[0,0,0]) %这里[0,0,0]为c的预估值，作为迭代初值。
MATLAB数学实验
第四章 函数和方程
第四章 函数和方程
4.1 预备知识：零点、极值和最小二乘法 4.2 函数零点、极值和最小二乘拟合的 MATLAB指令 4.3 计算实验：迭代法 4.4 建模实验：购房贷款的利率
4.1 预备知识：零点
非线性方程 f (x) = 0 • 若对于数有f () = 0, 则称为方程的解或根，也 称为函数f (x)的零点 • 若f () = 0, f ’()0 则称为单根。 • 若有k >1, f () = f ’() = …= f (k-1)() = 0，但 f (k)()0 , 称为k重根 • 非线性方程求解通常用数值方法求近似解 非线性方程（组）

男子110栏问题的求解
13.4 13.3 13.2 13.1 13 12.9 12.8 12.7 12.6 12.5 12.4 1970 1980 1990 2000 2010 2020 2030 2040 2050
大约2043年 突破12秒50
4.3 计算实验：迭代法
1 迭代法
迭代法是从解的初始近似值x0(简称初值)开 始，利用某种迭代格式x k+1 = g (x k ), 求得一近似值序列x1, x2, …, xk, xk+1, … 逐步逼近于所求的解(称为不动点)。 最常用的迭代法是牛顿迭代法，其迭代格式 为
2 ( y f ( c , x )) i i i 0 n
• 当f关于c是线性函数,问题转化为一个线性方程组求解， 且其解存在唯一。 • 如果f关于c是非线性函数，问题转化为函数极值问题
c= lsqnonlin (Fun,c0) 使用迭代法搜索最 优参数c. 其中Fun是以参数c(可以是向量) 为自变量的函数,表示误差向量 y-f(c,x)(x, y为数据)， c0为参数c的近似初值(与c同维向量)

男子110栏问题的求解
%考虑负指数函数，Matlab程序 clear;close all; data=[13.24 13.21 13.16 13.00 12.93 12.92 12.91 12.88 12.87]; year=[1972 1977 1979 1979 1981 1989 1993 2006 2008]; plot(year,data,'o'); y=data;x=year-1972; fun=@(c,x)c(1)*exp(-c(2)*x); c=lsqcurvefit(fun,[13,(13.24-12.87)/36],x,y) fun2=@(x)c(1)*exp(-c(2)*x)-12.5; x2=fsolve(fun2,50)%约71年以后，即2043年 xx=0:75; yy=fun(c,xx); figure; plot(year,data,'o',xx+1972,yy);grid on;
x k 1 x k
f ( xk ) f ' ( xk )
Newton法几何意义：切线法
切线代替曲线 0 f ( x ) f ( x
) f ( x )( x xk 1 ) k 1 k 1
f ( xk 1 ) xk xk 1 f ( xk 1 )
例3 求函数y=xsin(x2-x-1)在(-2, -0.1)内的 零点
>> fun=inline('x*sin(x^2-x-1)','x')
>> fzero(fun,[-2 -0.1])
>> fplot(fun,[-2,-0.1]),grid on;
>> fzero(fun,[-2,-1.2]), fzero(fun,[-1.2,-0.1]) >> fzero(fun,-1.6), fzero(fun,-0.6)
计算结果 c= 1.7427 -1.6958 1.0850 注意lsqnonlin与lsqcurvefit用法上的区别
最小二乘拟合

函数类型的选择
根据物理或数学机理 作图看趋势 常用直线、二次函数，指数函数等


百度文库
参数初始值的选择
根据参数的实际意义 由数据估计 网格化搜索 常用原点