MATLAB 程式設計進階篇 多項式的處理與分析

Presentation on theme: "MATLAB 程式設計進階篇 多項式的處理與分析"— Presentation transcript:

1 MATLAB 程式設計進階篇 多項式的處理與分析
張智星 清大資工系 多媒體檢索實驗室

2 7-1多項式的加減乘除 n 次多項式的表示法 MATLAB 中，用 n+1 的列向量表示 p(x)
Ex. p = [1,2,3,1] 表示一個三次多項式

3 多項式的加減 (I) 計算其和與差 多項式的加減，直接由向量的加減推出 兩多項式， 及 範例7-1：polyPlus.m
兩多項式， 及 計算其和與差 範例7-1：polyPlus.m p1 = [1,0,1,1]; p2 = [0,1,-1,2]; p1 + p2 p1 - p2

4 多項式的加減 (II) 矩陣 p1 與 p2 的長度要一致，否則 MATLAB 會產生運算錯誤的訊息 ans = 1 1 0 3
矩陣 p1 與 p2 的長度要一致，否則 MATLAB 會產生運算錯誤的訊息

5 多項式的乘除 (I) 多項式的乘除，使用 conv 及 deconv 指令達成 欲求多項式 與 的乘積 以多項式的表示法，為
欲求多項式 與 的乘積 範例7-2：conv01.m 以多項式的表示法，為 p1 = [1, 0, 1, 1]; p2 = [0, 1, -1, 2]; p3 = conv(p1, p2) p3 =

6 多項式的乘除 (II) 欲求 除以 的商式與餘式 代表 除以 ，得到商式 ，餘式 範例7-3：deconv01.m
欲求 除以 的商式與餘式 範例7-3：deconv01.m 代表 除以 ，得到商式 ，餘式 p1 = [1, 0, 1, 1]; p2 = [1, -1, 2]; [q, r] = deconv(p1, p2) q = r =

7 多項式的加減乘除 函數 說明 p1 + p2 多項式 與 的和 p1 - p2 多項式 與 的差 conv(p1, p2)
多項式 與 的和 p1 - p2 多項式 與 的差 conv(p1, p2) 多項式 與 的乘積 [q, r] = deconv(p1, p2) 多項式 除以 ，得到商式為 ，餘式為

8 7-2 多項式的求值、求根 、微分與積分 計算多項式的值，使用 polyval 指令 範例7-4：polyval01.m
x = 0:0.1:3; y = polyval(p, x); plot(x, y, '-o');

9 計算多項式的值-1 x 和 y 是長度為 31 的向量，y(i) 的值為　　　　　 在 x = x(i) 的函數值

10 計算多項式的值-2 計算p(A)，A為方陣，可用polyvalm指令 範例7-5：polyvalm01.m
結果和 B = A^2 + 2*A + 1 相同 p = [1 2 1]; A = [1 2; 3 4]; B = polyvalm(p, A) B =

11 計算多項式的根 (I) 多項式的根，可用 roots 指令求得 多項式 的根 範例7-6：roots01.m
多項式 的根 範例7-6：roots01.m p = [1, 3, 1, 5, -1]; % 多項式 r = roots(p); % 求多項式的根 r = i i 0.1886

12 計算多項式的根 (II) 為驗證此四根，可輸入如下 範例7-7：roots02.m p = [1, 3, 1, 5, -1]; % 多項式
r = roots(p); % 求多項式的根 polyval(p, r); % 將根帶入多項式求值 ans = 1.0e-012 * 0.1465 i i 0.0002

13 提示 fzero指令： Roots指令： 用於一般函數的求根 一次只找一個根 使用牛頓法 只用於多項式求根 一次找到全部的根
將多項式表示成「伴隨矩陣」（Companion Matrix），再用解特徵值方法求根

14 多項式的微分 求多項式的微分，使用指令polyder 範例7-8：polyder01.m 表示 微分後的結果為 p = [1 3 3 1];
表示 微分後的結果為 p = [ ]; q = polyder(p) q =

15 多項式的積分 (I) MATLAB 6.x 以後已經提供 polyint 指令 範例7-9：polyint02.m 表示 積分後為
表示 積分後為 p = [ ]; k = 8; % 積分後的不定常數 q = polyint(p, k) % 積分後的多項式 q =

16 多項式的積分 (II) MATLAB 5.x 並無對多項式積分的指令 用其它方法達成積分 範例7-10：polyint01.m
t = length(p):-1:1; k = 8; % 積分後的不定常數 q = [p./t, k] % 積分後的多項式 q =

17 多項式的值、根、積分、微分 函數 說明 q = polyval(p, x) 計算 p(x) 的值。 q = polyvalm(a, A)
計算 p(A)，A 為一方陣。 q = polyder(p) q(x) 為 p(x) 的微分。 q = [p./length(p):-1:1, k] q(x) 為 p(x) 的積分，其中 k 為任意常數。 q = polyint(p, k) q(x) 為 p(x) 的積分，其中 k 為任意常數。(本指令只適用於 MATLAB 6.x 以後的版本。)

18 矩陣的特徵多項式 給定一方陣 A，其特徵多項式為 使用 poly 指令計算特徵多項式 即為矩陣 A 的特徵多項式
範例7-11：poly01.m 即為矩陣 A 的特徵多項式 A = [1 3 4; 2 4 1; 1 6 2]; p = poly(A) p =

19 矩陣的特徵值 特徵多項式的根即為矩陣 A 的特徵值 範例7-12：poly02.m A = [1 3 4; 2 4 1; 1 6 2];
p = poly(A); % 方陣的特徵多項式 r = roots(p); % 特徵方程式的根，亦即固有值 det(A-r(1)*eye(3)) det(A-r(2)*eye(3)) det(A-r(3)*eye(3)) ans = e-014 ans = e e-015i ans = e e-015i

20 提示 特徵方程式的根會讓 ，也是方陣 A 的特徵值或固有值（Eigenvalues）的定義 可用eig指令來計算方陣 A 的特徵值及特徵向量

21 7-4 部份分式展開 指令 residue 計算一分式的部份分式展開 A(s) 和 B(s) 為多項式，且 B(s) 無重根
分式 可以表示為 , ,…, 為 的根(或 的極點) , ,…, 為常數 為一多項式

22 多項式部份分式展開 (I) 欲求 的部份分式展開 範例7-13：residue01.m b = [3 8]; a = [1 5 6];
欲求 的部份分式展開 範例7-13：residue01.m b = [3 8]; a = [1 5 6]; [r, p, k] = residue(b, a)

23 多項式部份分式展開 (II) 由以結果得知， 。 r = 1.0000 2.0000 p = k = []

24 提示 部份分式展開特別適用於線性系統之轉換函數（Transfer Function）的分析
若 為一系統之轉換函數，經由上述部份分式展開後，可知其脈衝響應（Impulse Response）為

25 多項式部份分式還原 residue 指令將部份分式轉回原來的形式 範例7-14：residue02.m 將部份分式結合成單一分式：
p = [-3 -2]; k = []; [b2, a2] = residue(r, p, k) b2 = a2 =

26 7-5多項式擬合 曲線擬合(Curve Fitting) 資料分析上常用的方法 以一參數化的曲線逼近一組給定的資料點
若參數化的曲線是一多項式，此曲線擬合稱為多項式擬合（Polynomial Fitting） polyfit 指令用以進行多項式擬合

27 多項式擬合範例 (I) 美國人口在1790至1990(10年一期)資料 範例7-15：plotCensus.m
load census.mat plot(cdate, pop, 'o'); xlabel('年份'); ylabel('人口（單位：百萬）');

28 多項式擬合範例 (II)

29 多項式擬合範例 (III) 欲預測美國在西元 2002 年的總人口，使用 polyfit 找出多項式逼近已知的資料點
範例7-16：polyfit01.m load census.mat plot(cdate, pop, 'o'); p3 = polyfit(cdate, pop, 3); cdate2 = 1790:2002; pop2 = polyval(p3, cdate2); plot(cdate, pop, 'o', cdate2, pop2, '-'); xlabel('年份'); ylabel('人口（單位：百萬）'); legend('實際人口', '預測人口'); popAt2002 = polyval(p3, 2002)

30 多項式擬合範例 (IV) 此曲線預測在 2002 年美國人口約為 百萬人 population =

31 不同次數的比較 (I) 不同的次數進行多項式擬合，「外插」（Extrapolation）時得到的差異相當大
範例7-17：polyfit02.m load census.mat cdate2 = min(cdate):(max(cdate)+30); curve = zeros(7, length(cdate2)); for i = 1:7 curve(i,:) = polyval(polyfit(cdate,pop,i), cdate2); end plot(cdate, pop,'o', cdate2, curve); legend('實際資料', 'order=1', 'order=2','order=3', 'order=4', 'order=5', 'order=6', 'order=7'); xlabel('年份'); ylabel('人口（單位：百萬）');

32 不同次數的比較 (II) 次數越高，越能逼近給定的資料 不表示預測的準確度會提高
高次多項式對雜訊（Noise）敏感度較高，容易產生不準確的外插預測

33 多項式擬合的誤差及錯誤 上述範例中，常出現 MATLAB 警告訊息 表示實際資料的不足 或資料數值範圍的過大或過小
在進行多項式擬合時，最好先將資料進行必要的平移、放大或縮小範圍