系統(tǒng)辨識—最小二乘法

1、最小二乘法參數(shù)辨識1 引言系統(tǒng)辨識是根據(jù)系統(tǒng)的輸入輸出時間函數(shù)來確定描述系統(tǒng)行為的數(shù)學(xué)模型?，F(xiàn)代控制理論中的一個分支。通過辨識建立數(shù)學(xué)模型的目的是估計表征系統(tǒng)行為的重要參數(shù)，建立一個能模仿真實系統(tǒng)行為的模型，用當(dāng)前可測量的系統(tǒng)的輸入和輸出預(yù)測系統(tǒng)輸出的未來演變，以及設(shè)計控制器。對系統(tǒng)進(jìn)行分析的主要問題是根據(jù)輸入時間函數(shù)和系統(tǒng)的特性來確定輸出信號。對系統(tǒng)進(jìn)行控制的主要問題是根據(jù)系統(tǒng)的特性設(shè)計控制輸入，使輸出滿足預(yù)先規(guī)定的要求。而系統(tǒng)辨識所研究的問題恰好是這些問題的逆問題。通常，預(yù)先給定一個模型類=M（即給定一類已知結(jié)構(gòu)的模型），一類輸入信號u和等價準(zhǔn)則J=L(y，yM)(一般情況下，J是誤差函數(shù)

2、，是過程輸出y和模型輸出yM的一個泛函)；然后選擇使誤差函數(shù)J達(dá)到最小的模型，作為辨識所要求的結(jié)果。系統(tǒng)辨識包括兩個方面：結(jié)構(gòu)辨識和參數(shù)估計。在實際的辨識過程中，隨著使用的方法不同，結(jié)構(gòu)辨識和參數(shù)估計這兩個方面并不是截然分開的，而是可以交織在一起進(jìn)行的。2 系統(tǒng)辨識的目的在提出和解決一個辨識問題時，明確最終使用模型的目的是至關(guān)重要的。它對模型類（模型結(jié)構(gòu)）、輸入信號和等價準(zhǔn)則的選擇都有很大的影響。通過辨識建立數(shù)學(xué)模型通常有四個目的。 估計具有特定物理意義的參數(shù)有些表征系統(tǒng)行為的重要參數(shù)是難以直接測量的，例如在生理、生態(tài)、環(huán)境、經(jīng)濟(jì)等系統(tǒng)中就常有這種情況。這就需要通過能觀測到的輸入輸出數(shù)據(jù)，用辨

3、識的方法去估計那些參數(shù)。 仿真仿真的核心是要建立一個能模仿真實系統(tǒng)行為的模型。用于系統(tǒng)分析的仿真模型要求能真實反映系統(tǒng)的特性。用于系統(tǒng)設(shè)計的仿真，則強調(diào)設(shè)計參數(shù)能正確地符合它本身的物理意義。 預(yù)測這是辨識的一個重要應(yīng)用方面，其目的是用迄今為止系統(tǒng)的可測量的輸入和輸出去預(yù)測系統(tǒng)輸出的未來的演變。例如最常見的氣象預(yù)報，洪水預(yù)報，其他如太陽黑子預(yù)報，市場價格的預(yù)測，河流污染物含量的預(yù)測等。預(yù)測模型辨識的等價準(zhǔn)則主要是使預(yù)測誤差平方和最小。只要預(yù)測誤差小就是好的預(yù)測模型，對模型的結(jié)構(gòu)及參數(shù)則很少再有其他要求。這時辨識的準(zhǔn)則和模型應(yīng)用的目的是一致的，因此可以得到較好的預(yù)測模型。 控制為了設(shè)計控制系統(tǒng)就需

4、要知道描述系統(tǒng)動態(tài)特性的數(shù)學(xué)模型，建立這些模型的目的在于設(shè)計控制器。建立什么樣的模型合適，取決于設(shè)計的方法和準(zhǔn)備采用的控制策略。3 系統(tǒng)辨識的方法經(jīng)典方法：經(jīng)典的系統(tǒng)辨識方法的發(fā)展已經(jīng)比較成熟和完善，他包括階躍響應(yīng)法、脈沖響應(yīng)法、頻率響應(yīng)法、相關(guān)分析法、譜分析法、最小二乘法和極大似然法等。其中最小二乘法(LS)是一種經(jīng)典的和最基本的，也是應(yīng)用最廣泛的方法。但是，最小二乘估計是非一致的，是有偏差的，所以為了克服他的缺陷，而形成了一些以最小二乘法為基礎(chǔ)的系統(tǒng)辨識方法：廣義最小二乘法(GI S)、輔助變量法(IV)、增廣最小二乘法(EI，S)和廣義最小二乘法(GI S)，以及將一般的最小二乘法與其他

5、方法相結(jié)合的方法，有最小二乘兩步法(CORI S)和隨機逼近算法等。經(jīng)典的系統(tǒng)辨識方法還存在著一定的不足： (1)利用最小二乘法的系統(tǒng)辨識法一般要求輸入信號已知，并且必須具有較豐富的變化，然而，這一點在某些動態(tài)系統(tǒng)中，系統(tǒng)的輸入常常無法保證；(2)極大似然法計算耗費大，可能得到的是損失函數(shù)的局部極小值；(3)經(jīng)典的辨識方法對于某些復(fù)雜系統(tǒng)在一些情況下無能為力?，F(xiàn)代方法：隨著系統(tǒng)的復(fù)雜化和對模型精確度要求的提高，系統(tǒng)辨識方法在不斷發(fā)展，特別是非線性系統(tǒng)辨識方法。主要有：1、集員系統(tǒng)辨識法在1979年集員辨識首先出現(xiàn)于Fogel 撰寫的文獻(xiàn)中，1982年Fogel和Huang又對其做了進(jìn)一步的改進(jìn)

6、。集員辨識是假設(shè)在噪聲或噪聲功率未知但有界UBB(Unknown But Bounded)的情況下，利用數(shù)據(jù)提供的信息給參數(shù)或傳遞函數(shù)確定一個總是包含真參數(shù)或傳遞函數(shù)的成員集(例如橢球體、多面體、平行六邊體等)。不同的實際應(yīng)用對象，集員成員集的定義也不同。集員辨識理論已廣泛應(yīng)用到多傳感器信息融合處理、軟測量技術(shù)、通訊、信號處理、魯棒控制及故障檢測等方面。2、多層遞階系統(tǒng)辨識法多層遞階方法的主要思想為：以時變參數(shù)模型的辨識方法作為基礎(chǔ)，在輸入輸出等價的意義下，把一大類非線性模型化為多層線性模型，為非線性系統(tǒng)的建模給出了一個十分有效的途徑。3、神經(jīng)網(wǎng)絡(luò)系統(tǒng)辨識法由于人工神經(jīng)網(wǎng)絡(luò)具有良好的非線性映射

7、能力、自學(xué)習(xí)適應(yīng)能力和并行信息處理能力，為解決未知不確定非線性系統(tǒng)的辨識問題提供了一條新的思路。與傳統(tǒng)的基于算法的辨識方法相比較，人工神經(jīng)網(wǎng)絡(luò)用于系統(tǒng)辨識具有以下優(yōu)點：(1)不要求建立實際系統(tǒng)的辨識格式，可以省去對系統(tǒng)建模這一步驟；(2)可以對本質(zhì)非線性系統(tǒng)進(jìn)行辨識；(3)辨識的收斂速度僅與神經(jīng)網(wǎng)絡(luò)的本身及所采用的學(xué)習(xí)算法有關(guān)；(4)通過調(diào)節(jié)神經(jīng)元之間的連接權(quán)即可使網(wǎng)絡(luò)的輸出來逼近系統(tǒng)的輸出；(5)神經(jīng)網(wǎng)絡(luò)也是系統(tǒng)的一個物理實現(xiàn)，可以用在在線控制。4、模糊邏輯系統(tǒng)辨識法模糊邏輯理論用模糊集合理論，從系統(tǒng)輸入和輸出的量測值來辨識系統(tǒng)的模糊模型，也是系統(tǒng)辨識的一個新的和有效的方法，在非線性系統(tǒng)辨識

8、領(lǐng)域中有十分廣泛的應(yīng)用。模糊邏輯辨識具有獨特的優(yōu)越性：能夠有效地辨識復(fù)雜和病態(tài)結(jié)構(gòu)的系統(tǒng)；能夠有效地辨識具有大時延、時變、多輸入單輸出的非線性復(fù)雜系統(tǒng)；可以辨識性能優(yōu)越的人類控制器；可以得到被控對象的定性與定量相結(jié)合的模型。模糊邏輯建模方法的主要內(nèi)容可分為兩個層次：一是模型結(jié)構(gòu)的辨識，另一個是模型參數(shù)的估計。典型的模糊結(jié)構(gòu)辨識方法有：模糊網(wǎng)格法、自適應(yīng)模糊網(wǎng)格法、模糊聚類法及模糊搜索樹法等。5、小波網(wǎng)絡(luò)系統(tǒng)辨識法小波網(wǎng)絡(luò)是在小波分解的基礎(chǔ)上提出的一種前饋神經(jīng)網(wǎng)絡(luò)口 ，使用小波網(wǎng)絡(luò)進(jìn)行動態(tài)系統(tǒng)辨識，成為神經(jīng)網(wǎng)絡(luò)辨識的一種新的方法。小波分析在理論上保證了小波網(wǎng)絡(luò)在非線性函數(shù)逼近中所具有的快速性、準(zhǔn)

9、確性和全局收斂性等優(yōu)點。小波理論在系統(tǒng)辨識中，尤其在非線性系統(tǒng)辨識中的應(yīng)用潛力越來越大，為不確定的復(fù)雜的非線性系統(tǒng)辨識提供了一種新的有效途徑，其具有良好的應(yīng)用前景。1. 設(shè)計內(nèi)容設(shè)SISO系統(tǒng)的差分方程為： 式（1-1）參數(shù)取真值為：，利用MATLAB 的M語言辨識系統(tǒng)中的未知參數(shù)、。2. 設(shè)計過程2.1 建立系統(tǒng)設(shè)SISO系統(tǒng)的差分方程為： 式（2-1）參數(shù)取真值為：，利用MATLAB 的M語言辨識系統(tǒng)中的未知參數(shù)、。要求：用參數(shù)的真值利用差分方程求出作為測量值，是均值為0，方差為0.1、0.5和0.01的不相關(guān)隨機序列。選取一種最小二乘算法辨識。2.2 最小二乘簡介2.2.1 最小二乘法的

10、概念與應(yīng)用對工程實踐中測得的數(shù)據(jù)進(jìn)行理論分析，用恰當(dāng)?shù)暮瘮?shù)去模擬數(shù)據(jù)原型是一類十分重要的問題，最常用的逼近原則是讓實測數(shù)據(jù)和估計數(shù)據(jù)之間的距離平方和最小，這即是最小二乘法。最小二乘法是一種經(jīng)典的數(shù)據(jù)處理方法。在系統(tǒng)辨識領(lǐng)域中 ,最小二乘法是一種得到廣泛應(yīng)用的估計方法 ,可用于動態(tài)系統(tǒng) ,靜態(tài)系統(tǒng) , 線性系統(tǒng) ,非線性系統(tǒng)?？捎糜陔x線估計，也可用于在線估計。這種辨識方法主要用于在線辨識。在隨機的環(huán)境下，利用最小二乘法時，并不要求觀測數(shù)據(jù)提供其概率統(tǒng)計方面的信息，而其估計結(jié)果，卻有相當(dāng)好的統(tǒng)計特性。MATLAB是一套高性能數(shù)字計算和可視化軟件 ,它集成概念設(shè)計 ,算法開發(fā) ,建模仿真 ,實時實現(xiàn)

11、于一體 ,構(gòu)成了一個使用方便、界面友好的用戶環(huán)境 ,其強大的擴(kuò)展功能為各領(lǐng)域的應(yīng)用提供了基礎(chǔ)。對于比較復(fù)雜的生產(chǎn)過程 ,由于過程的輸入輸出信號一般總是可以測量的 ,而且過程的動態(tài)特性必然表現(xiàn)在這些輸入輸出數(shù)據(jù)中 ,那么就可以利用輸入輸出數(shù)據(jù)所提供的信息來建立過程的數(shù)學(xué)模型。這種建模方法就稱為系統(tǒng)辨識。把辨識建模稱作“黑箱建模”。2.2.2 最小二乘法系統(tǒng)辨識結(jié)構(gòu)：本文把待辨識的過程看作“黑箱”。只考慮過程的輸入輸出特性，而不強調(diào)過程的內(nèi)部機理。+e(k) 圖1 SISO系統(tǒng)辨識“黑箱”結(jié)構(gòu)圖y(k)u(k)z(k)v(k)圖中，輸入u(k)和輸出z(k)是可以觀測的；G（）是系統(tǒng)模型，用來描述

12、系統(tǒng)的輸入輸出特性；N（）是噪聲模型，v(k)是白噪聲，e(k)是有色噪聲，根據(jù)表示定理：可以表示為e(k) =N（）v(k) 2.2.3 準(zhǔn)則函數(shù) 設(shè)一個隨機序列的均值是參數(shù)的線性函數(shù)： ，其中是可測的數(shù)據(jù)向量，那么利用隨機序列的一個實現(xiàn)，使準(zhǔn)則函數(shù)： (式2-2)達(dá)到極小的參數(shù)估計值稱作的最小二乘估計。最小二乘格式：，為模型參數(shù)向量，為零均值隨機噪聲。2.3 廣義最小二乘法2.3.1 廣義最小二乘數(shù)學(xué)模型式中，u(k)和表示系統(tǒng)的輸入輸出；v(k)是均值為零的不相關(guān)的隨機序列；且2.3.2 廣義最小二乘遞推算法如下式中2.3.3 廣義最小二乘遞推算法的計算步驟：1.給定初始條件 2利用式，

13、計算和；3利用式，構(gòu)造；4利用式遞推計算；5利用和計算；6根據(jù)來構(gòu)造；7利用返回第2步進(jìn)行迭代計算，直至獲得滿意的辨識結(jié)果。2.3.4 廣義最小二乘遞推算法的MATLAB仿真（程序源代碼見附錄）考慮仿真對象z(k)= -1.376z(k-1)-0.483z(k-2)+0.57u(k-1)+0.42u(k-2)+v(k)式中，v(k) 是均值為0，方差為0.01、0.1和0.5的不相關(guān)隨機序列。輸入信號采用4階M序列，幅度為1。選擇如下形式的辨識模型+y(k)u(k)e(k) z(k)+v(k)圖2廣義最小二乘法辨識實例結(jié)構(gòu)圖其中取c1=0,c2=0.3. 結(jié)果分析及算法優(yōu)化由于辨識算法中輸入或

14、噪聲信號為不相關(guān)隨機序列，所以每次辨識結(jié)果都不完全相同。但是，在相同輸入、相同的噪聲、相同的步長條件下，精度大體相同。算法優(yōu)化方案：（1）使用M序列（具有近似白噪聲的性質(zhì)）為輸入信號；（2）增加數(shù)據(jù)長度去L；（3）減小噪聲信號v(k)的方差。3.1 廣義最小二乘遞推算法的的MATLAB仿真結(jié)果及分析(1)、輸入選用題目給出的30個隨機數(shù)，即數(shù)據(jù)長度去L=30，噪聲選用均值為零，方差分別為0.5、0.1和0.01的隨機序列，辨識結(jié)果如表3-2-1表中給出了三種情況下辨識參數(shù)結(jié)果即表中的估計值，估計值與真值的相對誤差表3-2-1真值噪聲方差為0.5噪聲方差為0.1噪聲方差為0.01估計值相對誤差估

15、計值相對誤差估計值相對誤差a11.3761.5534-2.01051.47300.12311.37230.0063a20.4830.6884-2.0054 0.58920.1061 0.48150.0231b10.570.7916-2.02160.67900.1821 0.57990.0510b20.420.6213-2.0013 0.53610.1161 0.42180.0234輸入M序列，30步，噪聲方差0.5時： =12.556;輸入M序列，30步，噪聲方差0.1時： =2.5822;輸入M序列，30步，噪聲方差0.01時：=0.1706;（2）、輸入均采用M序列，噪聲選擇均值為零，方差

16、為0.5、0.1和0.01的隨機序列，辨識步長均為300步，辨識結(jié)果如表3-2-2。表中給出了三種情況下辨識參數(shù)結(jié)果即表中的估計值，估計值與真值的相對誤差. 表3-2-2真值噪聲方差為0.5噪聲方差為0.1噪聲方差為0.01估計值相對誤差估計值相對誤差估計值相對誤差a11.3761.3660-0.02801.37500.00791.37290.0005a20.4830.4749-0.0701 0.48590.0292 0.4710-0.0068b10.570.6613-0.12490.58270.0582 0.57200.0051b20.420.4112-0.0823 0.43440.0983

17、 0.4183-0.0049輸入采用M序列，噪聲方差0.5時： =135.288輸入采用M序列，噪聲方差0.1時： =28.917輸入采用M序列，噪聲方差0.01時：=2.7374（3）數(shù)據(jù)結(jié)果分析：輸入采用M序列比采用隨機序列得到的辨識效果更好。噪聲均值相等時，方差越大，辨識效果越差，反之，方差越小辨識效果越好。可以通過增加步長的方法提高辨識精度。 下面給出以M序列作為輸入，噪聲均值為零，方差為0.01的隨機序列，數(shù)據(jù)長度去L=30，得到的變化曲線圖： 下面給出以M序列作為輸入，噪聲均值為零，方差為0.01的不相關(guān)隨機序列，數(shù)據(jù)長度去L=300，得到的變化曲線圖：4. 參考文獻(xiàn)1) 侯媛彬,

18、 汪梅, 王立琦,系統(tǒng)辨識及其MATLAB仿真,科學(xué)出版社,2004年2) 方崇智，蕭德云，過程辨識，清華大學(xué)出版社，1988年3) 賈秋玲，袁冬莉，欒云鳳，基于MATLAB7.x/Simulink/Stateflow系統(tǒng)仿真、分析及設(shè)計，西北工業(yè)大學(xué)出版社，2006年4) 李言俊，張科，系統(tǒng)辨識理論及應(yīng)用，國防工業(yè)出版社，2006年課程心得上完這門課的第一次課，使我懂得了：辨識就是在輸入和輸出數(shù)據(jù)的基礎(chǔ)上，從一組給定的模型中，確定一個與所測系統(tǒng)等價的模型。了解了辨識的一些基本概念。知道了辨識的內(nèi)容和步驟。我們要判斷一個模型是否和原有系統(tǒng)等價或者近似，就好像我們要認(rèn)識一個人一樣，必須通過多種方

19、法，方式去測試他。這就需要我們給這模型充分的最優(yōu)激勵。這個形象的比喻讓我記憶深刻。接著我們學(xué)習(xí)了隨即過程的概念及相關(guān)知識，譜密度函數(shù)等內(nèi)容，這些為學(xué)習(xí)白噪聲，M序列等典型激勵提供了基礎(chǔ)。有了激勵U(k)必然有相應(yīng)的輸出Z(k)與之對應(yīng)，通過一批一批的數(shù)據(jù)U(k) 、Z(k)，通過用Matlab編程，能夠算出模型的參數(shù)。最后學(xué)習(xí)的是模型的辨識方法，辨識方法有很多，但根據(jù)模型形式可以分成兩類。一類是非參數(shù)模型辨識方法，另一種是參數(shù)模型辨識方法。第一種方法又稱經(jīng)典辨識方法，不必事先確定模型具體結(jié)構(gòu)，可適用于復(fù)雜的線性過程。第二種方法又稱現(xiàn)代辨識方法，必須假定一種結(jié)構(gòu)模型，通過極小化模型與過程之間的誤

20、差準(zhǔn)則函數(shù)來確定模型的參數(shù)。這種方法中我們主要學(xué)得是最小二乘類方法。這次作業(yè)選擇了一個較簡單的仿真實例，通過這次作業(yè)知道了一些函數(shù)和指令的使用方法。以后我想我會努力用一些時間來繼續(xù)學(xué)習(xí)這兩方面的知識內(nèi)容。最后感謝張老師的辛勤教育與諄諄教誨！附錄廣義最小二乘遞推算法的MATLAB仿真程序源代碼：clear %清理工作間變量L=300; % M序列的周期y1=1;y2=1;y3=1;y4=0; %四個移位寄存器的輸出初始值for i=1:L;%開始循環(huán)，長度為L x1=xor(y3,y4); %第一個移位寄存器的輸入是第三個與第四個移位寄存器的輸出的“或” x2=y1; %第二個移位寄存器的輸入是

21、第一個移位寄存器的輸出 x3=y2; %第三個移位寄存器的輸入是第二個移位寄存器的輸出 x4=y3; %第四個移位寄存器的輸入是第三個移位寄存器的輸出 y(i)=y4; %取出第四個移位寄存器的幅值為"0"和"1"的輸出信號，即M序列 if y(i)>0.5,u(i)=-1; %如果M序列的值為"1", 辨識的輸入信號取“-1” else u(i)=1; %如果M序列的值為"0", 辨識的輸入信號取“1” end %小循環(huán)結(jié)束 y1=x1;y2=x2;y3=x3;y4=x4; %為下一次的輸入信號做準(zhǔn)備end

22、 %大循環(huán)結(jié)束，產(chǎn)生輸入信號ufigure(1); %第一個圖形stem(u),grid on %顯示出輸入信號M序列徑線圖并給圖形加上網(wǎng)格v=normrnd(0, sqrt(0.01), 1, 300);%均值為零的，方差為0.01或0.5或0.1不相關(guān)的隨機噪聲ze(2)=0;ze(1)=0;for k=3:301; ze(k)=0*ze(k-1)+0*ze(k-2)+v(k-1);%C(z1)=1,即取c1=0,c2=0endz(2)=0;z(1)=0; %設(shè)z的前兩個初始值為零for k=3:301; %循環(huán)變量從3到301 z(k)=-1.376*z(k-1)-0.483*z(k-2

23、)+57*u(k-1)+0.42*u(k-2)+ze(k-1); %輸出采樣信號（測量值） end%RGLS廣義最小二乘辨識c0=0.0001 0.0001 0.0001 0.0001' %直接給出被辨識參數(shù)的初始值,即一個充分小的實向量pf0=106*eye(4,4); %直接給出初始狀態(tài)P0，即一個充分大的實數(shù)單位矩陣ce0=0.001 0.001'pe0=eye(2,2);c=c0,zeros(4,299); %被辨識參數(shù)矩陣的初始值及大小ce=ce0,zeros(2,299);e=zeros(4,300); %相對誤差的初始值及大小ee=zeros(2,300);s=0

24、;%廣義最小二乘遞推算法的計算步驟 for k=3:300; zf(k)=z(k)+ce(1,k-2)*z(k-1)+ce(2,k-2)*z(k-2); uf(k)=u(k)+ce(1,k-2)*u(k-1)+ce(2,k-2)*u(k-2); hf1=-zf(k-1),-zf(k-2),uf(k-1),uf(k-2)' x=hf1'*pf0*hf1+1; x1=inv(x); %開始求K(k) k1=pf0*hf1*x1;%求出K的值 d1=zf(k)-hf1'*c0; c1=c0+k1*d1; %求被辨識參數(shù)c e1=c1-c0; %求參數(shù)當(dāng)前值與上一次的值的差值 e2=e1./c0; %求參數(shù)的相對變化 e(:,k)=e2; %把當(dāng)前相對變化的列向量加入誤差矩陣的最后一列 c0=c1; %新獲得的參數(shù)作為下一次遞推的舊參數(shù) c(:,k)=c1; %把辨識參數(shù)c 列向量加入辨識參數(shù)矩陣的最后一列 pf1=pf0-k1*hf1'*pf0; %求出 p(k)的值 pf0=pf1; %給下次用 h1=-z(k-1),-z(k-2),u(k-1),u