Transcript 進階主題選讀

MATLAB 程式設計入門篇
進階主題選讀
修改自張智星教授講義
[email protected]
http://www.cs.nthu.edu.tw/~jang
清大資工系 多媒體資訊檢索實驗室
多項式加減乘除

n次多項式


在matlab中，用向量表示多項式



p(x) = anxn + an-1xn-1 + ... + a1x + a0
p = [an, an-1, ..., a1, a0]
例如[1 2 3 1]代表x3 + 2x2 + 3x + 1
運算：
p1 + p2
多項式和
p1 - p2
多項式差
conv(p1, p2)
多項式乘積
[q, r] = deconv(p1, p2) 多項式相除，得到商式和餘式
多項式求值、微積分

用polyval指令進行多項式求值：
p = [3 4 -2 5];
x = linspace(-3, 3);
y = polyval(p, x);
plot(x, y, '.-');


多項式微分：q = polyder(p)
多項式積分：q = polyint(p, k) % k為不定常數
多項式擬合

用polyfit指令進行多項式擬合：
% 製造資料
p = [2 -3 -4];
x = linspace(-2, 2, 20);
y = polyval(p, x) + (rand(1, 20)-0.5); % 加上適當誤差
% 擬合，假設已經知道是二次式
p2 = polyfit(x, y, 2);
y2 = polyval(p2, x);
% 畫圖
plot(x, y, 'o', x, y2, 'x-');
legend('原始資料', '擬合結果');
多項式擬合

擬合時的次方數差異：
% 製造資料
p = [2 -3 -4];
x = linspace(-1, 1, 20);
x2 = linspace(-1.3, 1.3);
y = polyval(p, x) + (rand(1, 20)-0.5); % 加上適當誤差
% 擬合，測試到七次方，擬合後之曲線有100個點，每個row是一條線
curve = zeros(7, 100);
for i = 1:7
thisCoef = polyfit(x, y, i);
curve(i, :) = polyval(thisCoef, x2);
end
% 畫圖
plot(x, y,'o', x2, curve);
legend('實際資料', '1次', '2次','3次', '4次', '5次', '6次', '7次');
單變數函數求根

範例：找出humps函數的根
 humps是一個單變數函數，MATLAB常用的測試函數之一
% 找某個點附近的根
x = fzero('humps', 1.5);
y = humps(x);
fprintf('x: %f, y:%f\n', x, y);
% 找某區間內的根。區間邊界的函數值必須異號，否則error
x = fzero('humps', [-1 1]);
y = humps(x);
fprintf('x: %f, y:%f\n', x, y);
% 顯示中間結果
opt = optimset('disp', 'iter');
x = fzero('humps', [-1, 1], opt);
fzero與roots的差異

fzero:




用於一般單變數函數的求根
一次只找一個根
使用牛頓法
roots:



只用於多項式求根
一次找到全部的根
將多項式表示成「伴隨矩陣」（Companion Matrix），再
用解特徵值方法求根
單變數函數最佳化

範例：找出humps函數的最小值
% 畫出圖形
fplot('humps', [-2 2]);
% 找某區間內的最小值
x = fminbnd('humps', 0, 1);
y = humps(x);
fprintf('x: %f, y:%f\n', x, y);
% 顯示中間結果
opt = optimset('disp', 'iter');
x = fminbnd('humps', 0, 1, opt);
% 放鬆誤差管制
opt = optimset('disp', 'iter', 'tolx', 0.1);
x = fminbnd('humps', 0, 1, opt);
多變數函數最佳化

範例：找出自訂函數的最小值

myFunc.m:
function y = myFunc(x)
y = 100*(x(2)-x(1)^2)^2+(1-x(1))^2;

fminsearch使用範例:
% 找某區間內的最小值
[x z]= fminsearch('myFunc', [-1.2, 1]);
fprintf('x: [%f %f], y:%f\n', x(1), x(2), z);
關於最佳化

optimset共有五十多個選項，其中





最佳化的結果，和起始點的選定有很大的關連。如果有好的
起始點



Display: 是否顯示中間過程
MaxIter: 最大疊代次數
TolX:終止搜尋的自變數容忍度
TolFun:終止搜尋的函數值容忍度
加快收斂速度
提高找到全域最佳解的機會
通常需要一再重複，才能收斂到最佳點
單變數函數的數值積分



quad: 適應式Simpson積分法(Adaptive Simpson Quadrature)
quadl: 適應式Lobatto積分法(Adaptive Lobatto Quadrature)
範例：
a = quad('x.^2 + 2*x + 1', -1, 1)
a = quad('exp(-x)', -2, 3)
一維內插: interp1


利用多項式進行內插運算
yi = interp1(x, y, x1, method)





x: 原始資料x座標
y: 原始資料y座標
xi: 要進行內插的x座標
method: 使用的方法
 'nearest': 鄰近點
 'linear': 線性
 'spline': 片段式三次spline
 'pchip' 或 'cubic': 保持型狀的片段式三次
xi必須落在x的範圍內
 否則應使用 yi = interp1(x, y, xi, method, 'extrap')
一維內插: interp1

範例：
x = 0:1:4*pi;
y = sin(x).*exp(-x/5);
xi = 0:0.1:4*pi;
y1 = interp1(x, y, xi, 'nearest');
y2 = interp1(x, y, xi, 'linear');
y3 = interp1(x, y, xi, 'pchip');
y4 = interp1(x, y, xi, 'spline');
plot(x, y, 'o', xi, y1, xi, y2, xi, y3, xi, y4);
legend('Original', 'Nearest', 'Linear', 'Pchip', 'Spline');