matlab在數(shù)值分析中的應(yīng)用學(xué)習(xí)教案

1、會計學(xué)1matlab在數(shù)值在數(shù)值(shz)分析中的應(yīng)用分析中的應(yīng)用第一頁，共41頁。27.1 多項式插值原理多項式插值原理7.2 分段分段(fn dun)線性插值線性插值7.3 離散數(shù)據(jù)的多項式擬合離散數(shù)據(jù)的多項式擬合7.3.1 線性最小二乘法擬合線性最小二乘法擬合7.3.2 非線性最小二乘法擬合非線性最小二乘法擬合第1頁/共41頁第二頁，共41頁。3 插值與擬合是來源于實際、又廣泛應(yīng)用于實際的兩種重要方法。隨著(su zhe)計算機的不斷發(fā)展及計算水平的不斷提高，它們已在國民生產(chǎn)和科學(xué)研究等方面扮演著越來越重要的角色。 第2頁/共41頁第三頁，共41頁。47.1 多項式插值原理多項式插值原理

2、(yunl)已知 f(x) 在點 xi 上的函數(shù)(hnsh)值 yi=f(xi), (i=0,1,2,n)求多項式 P(x)=a0 + a1x + anxn滿足(mnz): P(xk)= yk (k = 0,1,n)P10(t) f(t) f(x)第3頁/共41頁第四頁，共41頁。57.2 分段分段(fn dun)線性線性插值插值 所謂分段線性插值就是通過插值點用折線(zhxin)段連接起來逼近原曲線，這也是計算機繪制圖形的基本原理。實現(xiàn)分段線性插值不需編制函數(shù)程序，MATLAB自身提供了內(nèi)部函數(shù)interp1其主要用法如下：interp1(x,y,xi) 一維插值第4頁/共41頁第五頁，共4

3、1頁。6 yi=interp1(x,y,xi) 對一組點(x,y) 進行插值，計算插值點xi的函數(shù)(hnsh)值。x為節(jié)點向量值，y為對應(yīng)的節(jié)點函數(shù)(hnsh)值。如果y 為矩陣，則插值對y 的每一列進行，若y 的維數(shù)超出x 或 xi 的維數(shù)，則返回NaN。 yi=interp1(y,xi) 此格式默認(rèn)x=1:n ，n為向量y的元素個數(shù)值，或等于矩陣y的size(y,1)。第5頁/共41頁第六頁，共41頁。7 yi=interp1(x,y,xi,method) method用來指定插值的算法。默認(rèn)為線性算法。其值常用的可以是如下的字符串。 nearest 線性最近項插值。 linear 線性插

4、值。 spline 三次樣條插值。 cubic 三次插值。 pchip 分段(fn dun)三次Hermite插值。 所有的插值方法要求x是單調(diào)的。x 也可能并非連續(xù)等距的。第6頁/共41頁第七頁，共41頁。8正弦曲線的插值示例(shl)：x = 0:10; y = x.*sin(x); xx = 0:.25:10; yy = interp1(x,y,xx); plot(x,y,kd,xx,yy)第7頁/共41頁第八頁，共41頁。9例 產(chǎn)品產(chǎn)量擬合(n h)（離散數(shù)據(jù)的一維插值 ） year = 1900:10:2010; product = 75.995 91.972 105.711 123

5、.203 131.669 150.697 179.323 203.212 226.505 249.633 256.344 267.893 ;p1995 = interp1(year,product,1995)x = 1900:1:2010;y = interp1(year,product,x,pchip);plot(year,product,o,x,y)第8頁/共41頁第九頁，共41頁。10插值結(jié)果(ji gu)為：p1995 = 252.9885第9頁/共41頁第十頁，共41頁。11 Matlab也能夠完成二維插值的運算，相應(yīng)的函數(shù)為interp2，使用方法與interpl基本相同，只是輸入

6、和輸出(shch)的參數(shù)為矩陣，對應(yīng)于二維平面上的數(shù)據(jù)點。二維插值zi=interp2(x, y, z, xi, yi, method)第10頁/共41頁第十一頁，共41頁。12 返回矩陣ZI，其元素包含對應(yīng)于參量(cnling)XI與YI（可以是向量、或同型矩陣）的元素，即Zi(i,j)Xi(i,j),yi(i,j)。用戶可以輸入行向量和列向量Xi與Yi，此時，輸出向量Zi與矩陣meshgrid(xi,yi)是同型的。同時取決于由輸入矩陣X、Y與Z確定的二維函數(shù)Z=f(X,Y)。參量(cnling)X與Y必須是單調(diào)的，且相同的劃分格式，就像由命令meshgrid生成的一樣。若Xi與Yi中有在

7、X與Y范圍之外的點，則相應(yīng)地返回NAN（Not a Number）。第11頁/共41頁第十二頁，共41頁。13二維插值圖例：X,Y = meshgrid(-3:.25:3);Z = peaks(X,Y);XI,YI = meshgrid(-3:.125:3);ZZ = interp2(X,Y,Z,XI,YI);surfl(X,Y,Z);hold on;surfl(XI,YI,ZZ+15)axis(-3 3 -3 3 -5 20);shading flathold off第12頁/共41頁第十三頁，共41頁。14例years = 1950:10:1990;service = 10:10:30;w

8、age = 150.697 199.592 187.625 179.323 195.072 250.287 203.212 179.092 322.767 226.505 153.706 426.730 249.633 120.281 598.243;w = interp2(service,years,wage,15,1975)插值結(jié)果(ji gu)為：w =190.6288第13頁/共41頁第十四頁，共41頁。15三維插值三維插值vi = interp3(x,y,z,v,xi,yi,zi, method) 找出由參量找出由參量X,Y,Z決定的三元函數(shù)決定的三元函數(shù)V=V(X,Y,Z)在點（在

9、點（XI,YI,ZI）的值。參量）的值。參量XI,YI,ZI是同型陣列或向量。若向量參量是同型陣列或向量。若向量參量XI,YI,ZI是不同長度，不同方向（行或列）的是不同長度，不同方向（行或列）的向量，這時輸出參量向量，這時輸出參量VI與與Y1,Y2,Y3為同型矩陣為同型矩陣。其中。其中Y1,Y2,Y3為用命令為用命令meshgrid(XI,YI,ZI)生成的同型陣列。若插值點生成的同型陣列。若插值點(XI,YI,ZI)中有位中有位于點于點(X,Y,Z)之外的點，則相應(yīng)之外的點，則相應(yīng)(xingyng)地返地返回特殊變量值回特殊變量值NaN。第14頁/共41頁第十五頁，共41頁。16例x,y,

10、z,v = flow(20); xx,yy,zz = meshgrid(.1:.25:10, -3:.25:3, -3:.25:3);vv = interp3(x,y,z,v,xx,yy,zz); slice(xx,yy,zz,vv,6 9.5,1 2,-2 .2); shading interp;colormap cool 第15頁/共41頁第十六頁，共41頁。177.3 離散數(shù)據(jù)離散數(shù)據(jù)(shj)的多項式擬合的多項式擬合 在科學(xué)實驗(k xu sh yn)的統(tǒng)計方法研究中,往往要從一組實驗數(shù)據(jù)(xi, yi)中尋找出自變量x 和因變量y之間的函數(shù)關(guān)系y=f(x) 。 第16頁/共41頁第十

11、七頁，共41頁。18第17頁/共41頁第十八頁，共41頁。19第18頁/共41頁第十九頁，共41頁。20已知數(shù)據(jù)表 x x1 x2 xmf(x) y1 y2 ym求擬合函數(shù)(hnsh): (x) = a0 + a1x + + anxn使得mkknknkyxaxaa1210)(達(dá)到(d do)最小通俗通俗(tn s)的講：的講：第19頁/共41頁第二十頁，共41頁。217.3.1 線性最小二乘法線性最小二乘法(chngf)擬合擬合第20頁/共41頁第二十一頁，共41頁。22P=polyfit(X,Y,n) 求出(最小二乘意義(yy)下)n次擬合多項式P(x)=a0 xn + a1xn-1 + +

12、 an-1x + an計算結(jié)果為系數(shù): P= a0, a1, , an-1, an 在MATLAB中實現(xiàn)(shxin)最小二乘法擬合通常采用polyfit函數(shù)進行。 第21頁/共41頁第二十二頁，共41頁。23 函數(shù)polyfit是指用一個多項式函數(shù)來對已知數(shù)據(jù)進行擬合： x=0:0.1:1; y= -0.447 1.978 3.28 6.16 7.08 7.34 7.66 9.56 9.48 9.30 11.2 P= polyfit (x, y, 2)P= -9.8108 20.1293 -0.0317函數(shù)返回(fnhu)的是一個多項式系數(shù)的行向量，寫成多項式形式為： 0317. 01293

13、.208108. 92xx第22頁/共41頁第二十三頁，共41頁。24 為了使用polyfit，首先必須指定我們希望以多少階多項式對以上數(shù)據(jù)進行擬合，如果我們指定一階多項式，結(jié)果為線性近似，通常稱為線性回歸(hugu)。我們選擇二階多項式進行擬合。第23頁/共41頁第二十四頁，共41頁。25多項式求值命令(mng lng)：y1=polyval(P, x) 其中，P是n次多項式的系數(shù)，x是自變量的值，y1是多項式在x處的值。第24頁/共41頁第二十五頁，共41頁。26 為了比較擬合結(jié)果，我們繪制兩者的圖形： xi=linspace (0, 1, 100); %繪圖的X-軸數(shù)據(jù)。 Z=polyv

14、al (p, xi); %得到多項式在數(shù)據(jù)點處的值。 當(dāng)然(dngrn)，我們也可以選擇更高冪次的多項式進行擬合，如10階： p=polyfit (x, y, 10); xi=linspace (0, 1,100); z=ployval (p, xi); 第25頁/共41頁第二十六頁，共41頁。27 可以上機繪圖進行比較，曲線在數(shù)據(jù)點附近更加接近數(shù)據(jù)點的測量(cling)值了，但從整體上來說，曲線波動比較大，并不一定適合實際使用的需要，所以在進行高階曲線擬合時，“越高越好”的觀點不一定對的。 第26頁/共41頁第二十七頁，共41頁。28例例 血藥濃度問題血藥濃度問題 某人用快速靜脈注射某人用快

15、速靜脈注射(jn mi zh sh)方式一次注入藥物方式一次注入藥物300mg，在一定時刻采集，在一定時刻采集血樣，測得血藥濃度如下：血樣，測得血藥濃度如下：t=0.25,0.5,1,1.5,2,3,4,6,8c=19.21,18.15,15.36,14.10,12.89,9.32,7.45,5.24,3.01擬合血藥濃度隨時間的變化規(guī)律。擬合血藥濃度隨時間的變化規(guī)律。 第27頁/共41頁第二十八頁，共41頁。29第28頁/共41頁第二十九頁，共41頁。30用如下(rxi)Matlab命令（xueyao.m）t=0.25,0.5,1,1.5,2,3,4,6,8;c= 19.21,18.15,1

16、5.36,14.10,12.89,9.32,7.45,5.24,3.01;y=log(c);aa=polyfit(t,y,1);a=aa(1)b=aa(2)k=-ad=300;v=d/exp(b)cc=exp(b)*exp(a*t);plot(t,c,k+,t,cc,r) 第29頁/共41頁第三十頁，共41頁。31a = -0.2347b = 2.9943k = 0.2347v = 15.0219第30頁/共41頁第三十一頁，共41頁。327.3.2 非線性最小二乘法非線性最小二乘法(chngf)擬擬合合第31頁/共41頁第三十二頁，共41頁。33Matlab命令(mng lng)：（1）Ls

17、qnonlin(f,a0) 其中，f.m是描述函數(shù)R(a)的函數(shù)文件名，a0是初值。第32頁/共41頁第三十三頁，共41頁。34例例 血藥濃度問題血藥濃度問題(wnt)function f=ct(x)t=0.25,0.5,1,1.5,2,3,4,6,8;c= 19.21,18.15,15.36,14.10,12.89,9.32,7.45,5.24,3.01;f=c-x(1)*exp(x(2)*t);nongdu.mx0=10,0.5lsqnonlin(ct,x0) ans = 20.2413 -0.2420第33頁/共41頁第三十四頁，共41頁。35（2）lsqcurvefit(fun,a0,

18、x,y,LB,UB) 其中，f.m是描述(mio sh)函數(shù)f(x,a)的函數(shù)文件名，a0是初值，x=x1,xn，y=y1,yn為數(shù)據(jù)第34頁/共41頁第三十五頁，共41頁。36例例 血藥濃度問題血藥濃度問題(wnt)function f=xue(x,t)f=x(1)*exp(x(2)*t);xuenon.mt=0.25,0.5,1,1.5,2,3,4,6,8;c= 19.21,18.15,15.36,14.10,12.89,9.32,7.45,5.24,3.01;x0=10,0.5x=lsqcurvefit(xue,x0,t,c);yy=xue(x,t);plot(t,c,k*,t,yy,r) 第35頁/共41頁第三十六頁，共41頁。37第36頁/共41頁第三十七頁，共41頁。38例例 藥物的吸收與排除（口服）藥物的吸收與排除（口服） 假設(shè)口服劑量為假設(shè)口服劑量為d，中心室血藥濃度數(shù)據(jù)，中心室血藥濃度數(shù)據(jù)(shj)為為t=0.083,0.