數(shù)值分析第3章5-7節(jié).ppt

1 3 5曲線擬合的最小二乘法 3 5 1最小二乘法及其計(jì)算 在函數(shù)的最佳平方逼近中如果只在一組離散點(diǎn)集上給定 這就是科學(xué)實(shí)驗(yàn)中經(jīng)常見(jiàn)到的實(shí)驗(yàn)數(shù)據(jù)的曲線擬合 2 記誤差 則的各分量分別為個(gè)數(shù)據(jù)點(diǎn)上的誤差 3 設(shè)是上線性無(wú)關(guān)函數(shù)族 在中找一函數(shù) 使誤差平方和 5 1 這里 5 2 4 這個(gè)問(wèn)題稱為最小二乘逼近 幾何上稱為曲線擬合的最小二乘法 用最小二乘求擬合曲線時(shí) 首先要確定的形式 確定的形式問(wèn)題不僅是數(shù)學(xué)問(wèn)題 還與問(wèn)題的實(shí)際背景有關(guān) 通常要用問(wèn)題的運(yùn)動(dòng)規(guī)律及給定的數(shù)據(jù)進(jìn)行數(shù)據(jù)描圖 確定的形式 然后通過(guò)實(shí)際計(jì)算選出較好的結(jié)果 5 為了使問(wèn)題的提法更有一般性 通常在最小二乘法中考慮加權(quán)平方和 5 3 這里是上的權(quán)函數(shù) 它表示不同點(diǎn)處的數(shù)據(jù)比重不同 就是次多項(xiàng)式 若是次多項(xiàng)式 的一般表達(dá)式為 5 2 表示的線性形式 6 這樣 最小二乘問(wèn)題就轉(zhuǎn)化為求多元函數(shù) 5 4 的極小點(diǎn)問(wèn)題 由求多元函數(shù)極值的必要條件 有 使 5 3 取得最小 7 若記 5 5 上式可改寫(xiě)為 5 6 這個(gè)方程稱為法方程 可寫(xiě)成矩陣形式 8 其中 5 7 要使法方程 5 6 有惟一解 就要求矩陣非奇異 9 顯然在任意個(gè)點(diǎn)上滿足哈爾條件 函數(shù)的最小二乘解為 定義10 方程 5 6 存在惟一的解 從而得到 于是 10 這樣得到的 對(duì)任何形如 5 2 的 都有 故確是所求最小二乘解 11 一般可取 但這樣做當(dāng)時(shí) 通常對(duì)的簡(jiǎn)單情形都可通過(guò)求法方程 5 6 得到 給定的離散數(shù)據(jù) 例如 求解法方程 5 6 將出現(xiàn)系數(shù)矩陣為病態(tài)的問(wèn)題 若兩邊取對(duì)數(shù)得 12 例7 這樣就變成了形如 5 2 的線性模型 此時(shí) 若令 則 已知一組實(shí)驗(yàn)數(shù)據(jù)如下 求它的擬合曲線 13 解 從圖中看到各點(diǎn)在一條直線附近 故可選擇線性函數(shù)作擬合曲線 將所給數(shù)據(jù)在坐標(biāo)紙上標(biāo)出 見(jiàn)圖3 4 圖3 4 14 令 這里 故 15 解得 由 5 6 得方程組 于是所求擬合曲線為 16 關(guān)于多項(xiàng)式擬合 Matlab中有現(xiàn)成的程序 其中輸入?yún)?shù)為要擬合的數(shù)據(jù) 為擬合多項(xiàng)式的次數(shù) 輸出參數(shù)為擬合多項(xiàng)式的系數(shù) 利用下面的程序 可在Matlab中完成上例的多項(xiàng)式擬合 17 x 11233345 f 444 566688 5 aa poly x f 1 y polyval aa x plot x f r x y k xlabel x ylabel y gtext y s1 x 18 結(jié)果如下 19 例8 設(shè)數(shù)據(jù)由表3 1給出 用最小二乘法確定及 解 表中第4行為 通過(guò)描點(diǎn)可以看出數(shù)學(xué)模型為 它不是線性形式 用給定數(shù)據(jù)描圖可確定擬合曲線方程為 兩邊取對(duì)數(shù)得 20 若令 先將轉(zhuǎn)化為 為確定 根據(jù)最小二乘法 取 則得 數(shù)據(jù)表見(jiàn)表3 1 得 21 故有法方程 解得 于是得最小二乘擬合曲線為 22 利用下面的程序 可在Matlab中完成曲線擬合 x 1 001 251 501 752 00 y 5 105 796 537 458 46 y1 log y aa poly x y1 1 a aa 1 b exp aa 2 y2 b exp a x plot x y r x y2 k xlabel x ylabel y gtext y a exp bx 23 結(jié)果如下 24 3 5 2用正交多項(xiàng)式做最小二乘擬合 如果是關(guān)于點(diǎn)集 5 8 用最小二乘法得到的法方程組 5 6 其系數(shù)矩陣是病態(tài)的 25 5 9 則方程 5 6 的解為 且平方誤差為 26 接下來(lái)根據(jù)給定節(jié)點(diǎn)及權(quán)函數(shù) 構(gòu)造帶權(quán)正交的多項(xiàng)式 注意 用遞推公式表示 即 5 10 這里是首項(xiàng)系數(shù)為1的次多項(xiàng)式 27 5 11 下面用歸納法證明這樣給出的是正交的 28 假定對(duì)及 要證對(duì)均成立 由 5 10 有 由 5 10 第二式及 5 11 中的表達(dá)式 有 均成立 5 12 29 而 于是由 5 12 當(dāng)時(shí) 另外 是首項(xiàng)系數(shù)為1的次多項(xiàng)式 它可由 由歸納法假定 當(dāng)時(shí) 的線性組合表示 由歸納法假定又有 30 由假定有 再考慮 5 13 利用 5 11 中表達(dá)式及以上結(jié)果 得 31 至此已證明了由 5 10 及 5 11 確定的多項(xiàng)式 組成一個(gè)關(guān)于點(diǎn)集的正交系 用正交多項(xiàng)式的線性組合作最小二乘曲線擬合 只要根據(jù)公式 5 10 及 5 11 逐步求的同時(shí) 相應(yīng)計(jì)算出系數(shù) 最后 由和的表達(dá)式 5 11 有 32 并逐步把累加到中去 最后就可得到所求的 用這種方法編程序不用解方程組 只用遞推公式 并且當(dāng)逼近次數(shù)增加一次時(shí) 只要把程序中循環(huán)數(shù)加1 其余不用改變 這里可事先給定或在計(jì)算過(guò)程中根據(jù)誤差確定 擬合曲線 33 3 6最佳平方三角逼近與快速傅里葉變換 當(dāng)是周期函數(shù)時(shí) 顯然用三角多項(xiàng)式逼近比用代數(shù)多項(xiàng)式更合適 本節(jié)主要討論用三角多項(xiàng)式做最小平方逼近及快速傅里葉變換 簡(jiǎn)稱FFT算法 34 3 6 1最佳平方三角逼近與三角插值 設(shè)是以 為周期的平方可積函數(shù) 用三角多項(xiàng)式 6 1 作為最佳平方逼近函數(shù) 由于三角函數(shù)族 在上是正交函數(shù)族 于是在上的最小平方三角逼近多項(xiàng)式的系數(shù)是 35 稱為傅里葉系數(shù) 函數(shù)按傅里葉系數(shù)展開(kāi)得到的級(jí)數(shù) 6 3 就稱為傅里葉級(jí)數(shù) 6 2 36 只要在上分段連續(xù) 則級(jí)數(shù) 6 3 一致收斂到 對(duì)于最佳平方逼近多項(xiàng)式 6 1 有 由此可以得到相應(yīng)于 4 11 的貝塞爾不等式 因?yàn)橛疫叢灰蕾囉?左邊單調(diào)有界 所以級(jí)數(shù) 37 當(dāng)只在給定的離散點(diǎn)集 上已知時(shí) 則可類似得到離散點(diǎn)集正交性與相應(yīng)的離散傅里葉系數(shù) 下面只給出奇數(shù)個(gè)點(diǎn)的情形 收斂 并有 38 可以證明對(duì)任何成立 令 39 這表明函數(shù)族在點(diǎn)集 上正交 若令 其中 40 當(dāng)時(shí) 于是 6 4 就是三角插值多項(xiàng)式 系數(shù)仍由 6 4 表示 41 由于 所以函數(shù)族在區(qū)間上是正交的 42 當(dāng)時(shí) 個(gè)復(fù)向量具有如下正交性 6 5 43 事實(shí)上 令 于是 即 若 若 則有 則 從而 44 于是 若 這就證明了 6 5 成立 即是正交的 則 于是 因此 在個(gè)點(diǎn)上的最小二乘傅里葉逼近為 45 6 6 其中 6 7 于是由 6 6 得 即 6 8 46 而 6 8 是由求的過(guò)程 稱為反變換 47 3 6 2快速傅氏變換 FFT 不論是按 6 7 式由求 由求 6 9 其中 正變換 或 反變換 還是由 6 4 計(jì)算傅里葉逼近系數(shù) 都可歸結(jié)為計(jì)算 是已知復(fù)數(shù)序列 或是按 6 8 48 當(dāng)較大且處理數(shù)據(jù)很多時(shí) 就是用高速的電子計(jì)算機(jī) 很多實(shí)際問(wèn)題仍然無(wú)法計(jì)算 直到20世紀(jì)60年代中期產(chǎn)生了FFT算法 大大提高了運(yùn)算速度 從而使傅氏變換得以廣泛應(yīng)用 FFT算法的基本思想就是盡量減少乘法次數(shù) 49 用 6 9 計(jì)算全部 表面看要做個(gè)乘法 特別當(dāng)時(shí) 只有個(gè)不同的值 因此可把同一個(gè)對(duì)應(yīng)的相加后再乘 這就能大量減少乘法次數(shù) 50 設(shè)正整數(shù)除以后得商及余數(shù) 則 稱為的同余數(shù) 以表示 由于 因此計(jì)算時(shí)可用的同余數(shù)代替 從而推出FFT算法 以為例 說(shuō)明FFT的計(jì)算方法 由于則 6 9 的和是 6 10 故有 51 將用二進(jìn)制表示為 其中只能取0或1 例如 根據(jù)表示法 有 公式 6 10 可表示為 52 6 11 若引入記號(hào) 6 12 53 則 6 11 變成 它說(shuō)明利用同余數(shù)可把計(jì)算分為步 用公式 6 12 計(jì)算 每計(jì)算一個(gè)只用2次復(fù)數(shù)乘法 計(jì)算一個(gè)用 次復(fù)數(shù)乘法 計(jì)算全部共用次復(fù)數(shù)乘法 若注意公式 6 12 還可進(jìn)一步簡(jiǎn)化為 54 將這表達(dá)式中二進(jìn)制表示還原為十進(jìn)制表示 55 6 13 同樣 6 12 中的也可簡(jiǎn)化為 即 即得 56 把二進(jìn)制表示還原為十進(jìn)制表示 得 6 14 同理 6 12 中可簡(jiǎn)化為 即 57 表示為十進(jìn)制 有 6 15 58 根據(jù)公式 6 13 6 14 6 15 由 逐次計(jì)算到 見(jiàn)表3 2 略 上面推導(dǎo)的的計(jì)算公式可類似地推廣到的情形 根據(jù)公式 6 13 6 14 6 15 一般情況的FFT計(jì)算公式如下 59 6 16 其中 從出發(fā) 由到算到 一組占用個(gè)復(fù)數(shù)單元 計(jì)算時(shí)需給出兩組單元 括號(hào)內(nèi)的數(shù)代表它的位置 在計(jì)算機(jī)中代表存放數(shù)的地址 即為所求 60 這個(gè)計(jì)算公式除了具有不倒地址的優(yōu)點(diǎn)外 計(jì)算只有兩重循環(huán) 計(jì)算過(guò)程中只要按地址號(hào)存放則最后得到的 就是所求離散頻譜的次序 外循環(huán)由計(jì)算到 內(nèi)循環(huán)由計(jì)算到 由計(jì)算到 更重要的是整個(gè)計(jì)算過(guò)程省計(jì)算量 由公式看到算一個(gè)共做次復(fù)數(shù)乘法 而最后一步計(jì)算時(shí) 由于 61 當(dāng)時(shí)比值是它比一般FFT的計(jì)算量 次乘法 也快一倍 我們稱 6 16 的計(jì)算公式為改進(jìn)的FFT算法 62 3 7有理逼近 3 7 1有理逼近與連分式 有理函數(shù)逼近是指用形如 的函數(shù)逼近 與前面討論一樣 如果最小就可得到最佳有理一致逼近 7 1 63 如果最小則可得到最佳有理平方逼近函數(shù) 本節(jié)主要討論利用函數(shù)的泰勒展開(kāi)獲得有理逼近函數(shù)的方法 對(duì)函數(shù)用泰勒展開(kāi)得 7 2 取部分和 64 7 3 65 7 4 7 3 右端為的無(wú)窮連分式的前5項(xiàng) 最后式子 若取 7 3 的前2 4 6 8項(xiàng) 則可分別得到的以下有理逼近 是它的緊湊形式 66 若用同樣多項(xiàng)的泰勒展開(kāi)部分和逼近 并計(jì)算處的值及 計(jì)算結(jié)果見(jiàn)表3 2 的準(zhǔn)確值為 從表3 1可以看出 67 但它們的計(jì)算量是相當(dāng)?shù)?這說(shuō)明用有理逼近比多項(xiàng)式逼近好得多 由此看出的精度比高出近10萬(wàn)倍 例9 用輾轉(zhuǎn)相除法將它化為連分式并寫(xiě)成緊湊形式 解 給出有理函數(shù) 用輾轉(zhuǎn)相除可逐步得到 68 本例中用連分式計(jì)算的值只需3次除法 1次乘法和7次加法 69 若直接用多項(xiàng)式計(jì)算的秦九韶算法則需6次乘法和1次除法及7次加法 可見(jiàn)將化成連分式可節(jié)省計(jì)算乘除法次數(shù) 對(duì)一般的有理函數(shù) 7 1 可轉(zhuǎn)化為一個(gè)連分式 它的乘除法運(yùn)算只需次 而直接用有理函數(shù) 7 1 計(jì)算乘除法次數(shù)為次 70 3 7 2帕德逼近 利用函數(shù)的泰勒展開(kāi)可以得到它的有理逼近 設(shè)在的泰勒展開(kāi)為 7 5 它的部分和記作 7 6 71 定義11 設(shè) 其中無(wú)公因式 且滿足條件 7 8 則稱為函數(shù)在處的階帕德逼近 記作 簡(jiǎn)稱的帕德逼近 如果有理函數(shù) 7 7 72 根據(jù)定義 若令 則滿足條件 7 8 等價(jià)于 即 由于應(yīng)用萊布尼茨求導(dǎo)公式得 73 這里是由 7 6 得到的 上式兩端除 并由可得 7 9 及 7 10 注意當(dāng)時(shí) 故 7 10 可寫(xiě)成 74 7 11 其中時(shí) 若記 7 12 75 則方程組 7 11 的矩陣形式為 定理10 76 根據(jù)定理10 求的帕