现代信号处理课件

第一章離散隨機信號統計分析基礎本章目的：對隨機信號做一個簡短的回顧，並且介紹一些在以後各章中我們將要用的概念。主要包括：1。隨機信號的基本概念

2。平穩隨機信號的時域統計描述

3。平穩隨機信號的z域及頻域的統計描述

4。線性系統對隨機信號的回應

5。隨機信號的模型

隨機信號分析二、隨機信號的統計特性

要完整地描述一個各態歷經隨機過程，理論上要有無限長時間記錄。但實際上這是不可能的。通常用統計方法對以下三個方面進行數學描述：

1）幅值域描述:均值、方均值、方差、概率密度函數等。（2）時間域描述:自相關函數、互相關函數。（3）頻率域描述:自功率譜密度函數、互功率譜密度函數。自相關函數性質自相關函數的應用當延時

很大時，隨機雜訊的自相關函數趨於零，而週期信號的自相關函數仍是週期函數，且其週期不變。互相關函數描述一個信號的取值對另一個信號的依賴程度。互相關函數具有以下性質：

①兩週期信號具有相同的頻率，才有互相關函數，即兩個非同頻的週期信號是不相關的。

②兩個相同週期的信號的互相關函數仍是週期函數，其週期與原信號的週期相同，並不丟失相位資訊。

③兩信號錯開一個時間間隔

0處相關程度有可能最高，它反映兩信號x(t)、y（t)之間主傳輸通道的滯後時間。自相關函數和自譜密度函數構成一對傅立葉變換對。自譜密度函數是從頻率域對隨機過程作統計描述，集中顯示了隨機過程的頻率結構。實際應用中，-f不可能出現，所以往往處理成單邊譜。雙邊譜Sx(f)與單邊Gx(f)的關係為：自功率譜密度函數1-1隨機信號的基本概念1。隨機變數

:

離散的，連續的對隨機變數X，我們一般用它的分佈函數、概率密度等特徵來描述。概率分佈函數： 概率密度：

均值：

均方值：

方差： 常見的隨機變數概率密度均勻分佈：高斯分佈：又叫正態分佈，4個特點及圖2。隨機信號（以上有關隨機變數的描速方法可推廣到隨機信號）定義：圖1電晶體直流放大器的溫漂電壓（連續的）

如果我們把對溫漂電壓的觀察看作為一個隨機試驗，那麼，每一次的記錄，就是隨機試驗的一次實現，相應的結果就是一個樣本函數：

所有样本函数的集合i＝1,2,…,N,N→∞，就構成了溫漂電壓可能經歷的整個過程，該集合就是一個隨機過程，也即隨機信號，記之為：X(t)

物理意义：是一個隨機變數當t在時間軸上取值時，我們可得到m個隨機變數，顯然，描述這m個隨機變數最全面的方法是利用其m維的概率分佈函數（或概率密度）：m維的概率分佈函數

（理论上有意义，实际应用困难繁琐，一阶二阶）離散隨機信號：X(n)數學特徵是時間n的函數（物理意義）：1：均值（數學期望）2：方差3：均方值å=¥®==NiNXinxNnXEnm1),(1lim)]([)(（4）自相關函數 （5）自協方差函數：“集合平均”，該集合平均是由X(n)的無窮樣本在相應時刻對應相加（或相乘後再相加）來實現。

（6）互相關函數 （7）互協方差函數注：如果稱信號X和Y是不相關的。可得：

§1.2平穩隨機信號的時域統計描述1。定義：一個離散隨機信號X(n)，如果其均值與時間n無關，其自相關函數和的選取無關，而僅和之差有關，那麼，我們稱X(n)為寬平穩的隨機信號，或廣義平穩隨機信號。其具有以下的統計特徵.1）均值為常值。

2）自相關函數和自協方差函數均只是m的函數。目的：使問題簡化，實際工程中大部分屬於這種嚴平穩隨機信號：指概率特性不隨時間的平移而變化（或說與時間基準點無關）的隨機信號。只有當X(n)是高斯隨機過程時，寬平穩才是嚴平穩。2。平穩隨機信號相關函數的性質

性質1：

性質2：若X(n)是實信號，則，即自相關函數為實偶函數；若X(n)是復信號，則。即自相關函數是Hermitian對稱的。性質3

：性質4：

自相關函數和自協方差函數的關係

12當時工程實際中，當m趨於無窮大時，可以認為不相關，存在：自相關函數反映的其他資訊例1.2.1

隨機相位正弦序列式中A,f均為常數，Φ是一隨機變數，在0~2π內服從均勻分佈，即顯然，對應Φ的一個取值，可得到一條正弦曲線（因為Φ在0~2π內的取值是隨機的，所以其每一個樣本x(n)都是一條正弦信號）。求其均值及其自相關函數，並判斷其平穩性。

解由定義，X(n)的均值和自相關分別是：由於 及 所以隨機相位正弦波是寬平穩的。例1.2.2

隨機振幅正弦序列如下式所示：式f中為常數，A為正態隨機變數，A：N（0，σ2），試求X(n)的均值、自相關函數，並討論其平穩性。解:均值 對於給定的時刻n，為一常數，所以自相關函數 由此可以看出，雖然X(n)的均值和時間無關，但其自相關函數不能寫成的形式，也即和的選取位置有關，所以隨機振幅正弦波不是寬平穩的。3平穩隨機信號的各態遍曆性為什麼？為了要精確地求出，需要知道的無窮多個樣本，即，這在實際工作中顯然是不現實的。因為我們在實際工作中能得到的往往是對的一次實驗記錄，也即一個樣本函數。什麼是？定義、好處

對一平穩隨機信號，如果它的所有樣本函數在某一固定時刻的一階和二階統計特性和單一樣本函數在長時間內的統計特性一致，我們則稱為各態遍曆信號。應用：各態遍曆平穩隨機信號有：例1.2.3

討論例1.2.1隨機相位正弦序列的各態遍曆性。

解對，其單一的時間樣本，為一常數，對

作時間平均，顯然由於上式是n對求和，故求和號中的第一項與n無關，而第二項應等於零，所以這和例1.2.1按集合平均求出的結果一樣，所以隨機相位正弦波既是平穩的，也是各態遍曆的。例1.2.4

隨機信號的取值在（-1，1）之間均勻分佈，但對每一個樣本，其值不隨時間變化，如圖1.2.2所示，試討論其平穩性和各態遍曆性。解如圖所示，顯然的集合均值始終等於零，集合自相關也和的選取位置無關，因此它是寬平穩的。但對單一的樣本，它的時間均值並不等於零，因此，不是各態遍曆的。圖1.2.2例1.2.4中的X(n)

§1.3平穩隨機信號的z域及頻域的統計表達1.3.1與的z變換及其收斂域

＝=

1.3.2功率譜

令

有

=

=

=所以 = (1.3.10)因而有 =

=所以 = (1.3.10)因而有 =

=

=

=

白雜訊例：設一隨機相位余弦信號式中，A、ω。為常數，θ是在（0，2π）上均勻分佈的隨機變數，其概率密度為：解：①由自相關函數定義，有（週期函數的自相關函數也是週期函數，且具有相同的週期。）②由功率譜定義，得（其功率譜為兩個沖激函數）例：設隨機相位余弦信號§1.4線性系統對隨機信號的回應

要解决的问题：

已知输入信号的统计特性，其经过线性移不变系统后输出的统计特性

1.均值：（物理意義）1自相關函數及功率譜（假設非平穩的，證明是平穩的）令l=r-k，式（1.4.3）可表示為

結論：確定信號的輸出等於輸入與系統的單位沖激回應的卷積；

隨機信號的輸出（自相關函數）等於輸入信號（自相關函數）與系統的單位沖激回應（自相關函數）的卷積；也可描述為：卷積的相關，等於相關的卷積，稱為：相關－卷積定理

即可用公式描述為：

維納--辛欽定理

物理意義：一個隨機信號通過系統，從頻域看其輸出功率譜密度等於輸入功率譜密度與系統頻率回應的模平方的乘積。圖1.4.1理想帶通濾波器的頻率回應3。互相關函數和互功率譜密度例1.4.1

在圖1.4.2中，如果已知隨機信號x與y的互相關函數，試證明：（1）（2）（3）證明:（1）按定義所以

(1.4.21)又因為，所以式（1.4.21）將上二式代入得上式將z用z-1代入得(1.4.22)（2）同理可證所以

又因為，所以，式（1.4.22）成為上式將z用z-1代入得重要結論:假設d(k)和x(k)都是平穩隨機過程，且

是已知的

§1.5隨機信號的模型1.5.1ARMA模型1.5.2MA模型1.5.3AR模型1.5.4三種模型係數間關係本章重點1。隨機信號的統計特徵及物理意義2。平穩隨機信號的定義（特徵）3。各態歷經（遍曆）平穩隨機信號的特徵4。自（互）相關函數的定義及物理意義及與其它特徵量之間的關係。5。白雜訊的定義及特性應用6。隨機信號作用於線性系統的回應（時域、頻域（Z、W））：相關－卷積定理；維納--辛欽定理；在估計線性時不變系統的沖激回應或頻率回應中的應用。7。隨機信號的模型：AR、MA、ARMA間轉換的關係。共128頁82第二章維納濾波與卡爾曼濾波2.1引言

本章內容是解決雜訊中提取信號的問題。而維納(Wiener)濾波與卡爾曼(Kalman)濾波就是用來解決這樣一類從雜訊中提取信號問題的一種過濾(或濾波)方法。信號的恢復及信號的判決是兩個概念。而信號恢復（估計）的準確性判斷準則也不一樣。如：誤差的代數和最小；誤差的絕對值和最小；誤差的平方和最小等。（見實例）共128頁83共128頁84Wiener&KalmanFilter設計（即獲得系統的單位沖激相應）的準則：（條件）滿足最小均方誤差（正交性原理）為準則的，即保證：好處：運算簡單，對過大雜訊敏感，小雜訊不敏感；符合實際工程需要因此維納過濾與卡爾曼過濾又常常被稱為最佳線性過濾與預測或線性最優估計。共128頁85輸入信號及濾波器輸出模型輸入信號：輸出信號：共128頁86

維納濾波器的輸入—輸出關係實質是設計系統(傳遞)函數或單位沖激回應

共128頁872.2維納濾波器的離散形式(I)——

時域解h(n)=0當n<0，因果系統共128頁88共128頁89共128頁90

用幾何圖形理解正交性原理共128頁91上式展開有：共128頁92共128頁93因此維納濾波器的設計問題，歸結為求解wiener-Hopf方程共128頁94

解Wiener-Hopf方程

存在的問題及解決的思路存在的問題：1。假設設計的是因果系統，由於存在k>0的約束條件，卷積定理（雙邊Z變換）不能用。2。實際物理系統為因果系統。解決思路：1。設計一個非因果性系統（濾波器）。2。用有限長的因果序列h(n)來逼近hopt(n).實質為設計FIR型濾波器。共128頁95共128頁96共128頁97共128頁98共128頁99共128頁100

結論1.維納濾波器的設計實質為求解Weiner-Hopf方程。2。非因果系統設計簡單3。最小方差準則的維納濾波器，用有限沖激回應的FIR濾波器來實現，計算複雜，工作量大，並不是有效的辦法。共128頁101例題1：共128頁1022.3維納濾波器的z域解求解的基本思想：

把x(n)加以白化來求維納-霍夫方程的z域解.(這種方法是由波德(Bode)和香農(Shannon)首先提出的)

白化：任何具有有理功率譜密度的隨機信號都可以看成是由一白色雜訊激勵一物理網路所形成。

共128頁103

共128頁104圖2.3s(n)的信號模型圖2.4x(n)的信號模型圖2.5信號模型共128頁105

B(z)是由單位圓內的零極點組成，B(1/z)是由對應的單位圓外的零極點組成。因此：B（Z）是一個因果(或物理可實現)的並且是最小相移的網路共128頁106共128頁107

用白化x(n)的方法求解維納-霍夫方程共128頁108於是，求在最小均方誤差下的最佳Hopt(z)的問題就歸結為求最佳G(z)的問題了。我們可以對G(z)加以因果性或非因果性的約束具體求解。由於G(z)的激勵源是將x(n)白化後得到的白雜訊，這就使得求圖2.7(b)中的最佳G(z)比求圖2.7(a)中的最佳H(z)容易，（為什麼白化的原因）得：共128頁1092.3.1沒有物理可實現性約束的(非因果的)維納濾波器）共128頁110共128頁111共128頁112共128頁113共128頁114共128頁115共128頁116共128頁117物理意義及其解釋共128頁118圖2.8決定於與特性的例子共128頁119最小均方差求解共128頁120共128頁121共128頁122最小均方誤差不僅與輸入信號的功率有關（反），而且與信號和雜訊的功率譜的乘積有關（正）。共128頁123

2.3.2有物理可實現性約束的(因果的)維納濾波器共128頁124共128頁125共128頁126

[例1]

設已知，以及

(白雜訊)

其中s(n)代表所希望得到的信號，代表加性白雜訊。求物理可實現與物理不可實現這二種情況下的及。共128頁127

[解]

因為所以

又因為共128頁128其中B(z)由單位圓內的零極點組成，B(z^-1)由單位圓外的零極點組成，上兩式比較得

(1)物理可實現情況共128頁129因為對於項。所以共128頁130利用式(2.52)並考慮到，得共128頁131取單位圓為積分圍線，上式等於單位圓內的積點的留數之和，即共128頁132而在經過此濾波器以前的均方誤差為所以通過維納濾波器後均方誤差下降8/3(≈2.7)倍。共128頁133

(2)非物理可實現的情況共128頁134取單位圓為積分圍線C。在單位圓內有二個極點：共128頁135。H式等於該二個極點的留數，因此前面求得物理可實現的所以在此例中非物理可實現情況的均方誤差略小於(即稍好於)物理可實現的情況。可以證明，物理可實現情況的最小均方誤差總不會小於非物理可實現的情況。共128頁1362.4維納預測器

2.4.1預測的可能性共128頁137共128頁138

2.4.2預測器的計算公式維納預測器共128頁139共128頁140共128頁141

(1)沒有物理可實現性約束的(非因果的)維納預測器共128頁142共128頁143(2)有物理可實現性約束的(因果的)維納預測器共128頁144共128頁145

2.4.3純預測器(N步)共128頁146共128頁147共128頁148共128頁149

[例2]

已知及

求(1)使均方誤差最小的

(2)最小均方誤差 共128頁150

[解]

因為所以由式(2.63)得因果的維納預測器，應有共128頁151因為所以共128頁152對只取上式中n>0的部分，得再回到z-域，得共128頁153代入Hopt(z)運算式，得這個結果可用方框圖表示在圖2.11中。圖2.11純預測的例子共128頁154

由式(2.65)得它說明：N越大，誤差越大，如果N=0則沒有誤差。共128頁155

現在來解釋上述疑問。我們是把看成由白雜訊通過B(z)產生的，而故該信號模型可以用一個一階差分方程來表達：共128頁156圖2.12一階信號模型的例子共128頁157共128頁158共128頁159

2.4.4維納預測器的時域解——

一步線性預測公式共128頁160共128頁161共128頁162共128頁163共128頁164共128頁165共128頁166

§2.5卡爾曼濾波的信號模型——狀態方程與量測方程

2.5.1離散狀態方程及其解共128頁167共128頁168共128頁169

2.5.2量測方程共128頁170圖2.13卡爾曼濾波的信號模型(a)多維情況；(b)一維情況。共128頁171

[例3]

仍沿用前面維納濾波中的例子(例1)，設，已知求卡爾曼信號模型中的與共128頁172[解]因為

所以共128頁173變換到時域： 所以

又，因為，所以。共128頁174

§2.6卡爾曼過濾的方法與公式共128頁175共128頁176共128頁177共128頁178共128頁179共128頁180共128頁181共128頁182共128頁183共128頁184共128頁185共128頁186共128頁187共128頁188共128頁189共128頁190

[例4]

設為實離散時間隨機過程，具有功率譜密度：並已知在k=0時開始觀察信號，試用卡爾曼過濾的計算公式求，並將計算結果與維納過濾方法計算(例1)的結果進行比較(此例與例1中的相同)。共128頁191[解]從給定的，可以求得的狀態方程。因為(參見§2.5，例3)所以

共128頁192又由於所以

代入式(2.89)、(2.102)、(2.99)及(2.103)，得共128頁193(2.105)(2.106)(2.107)(2.108)將式(2.106)代入式(2.108)，得共128頁194

從式(2.107)與式(2.109)中消去，得求穩態解，將式(2.110)中的代入並化簡，得共128頁195所以(只取正值解)所以共128頁196所以由此可見，已知前一個估計值

與當前的量測值，就可以求得當前的估計值。

共128頁197共128頁198共128頁199共128頁200共128頁201共128頁202共128頁203本章小結1。掌握Wiener濾波器的作用－－雜訊中信號的恢復（估計）。2。掌握Wiener濾波器的設計原則－－保證估計結果的均方誤差最小（正交原理）。3。掌握Wiener濾波器設計實質及結果－－實質為解Wiener－Hofp方程，結果形式為單位沖激回應h(n)(時域)或系統函數H(z)（頻域）.共128頁204共128頁205

z=0.01*randn(1,150)-0.37727;xx(1)=-0.31;Q=0.000001;R=1;p(1)=0.02;s(1)=-0.37727;fork=2:1:150s(k)=-0.37727;xs(k)=xx(k-1);ps(k)=p(k-1)+Q;K(k)=ps(k)/(ps(k)+R);xx(k)=xs(k)+K(k)*(z(k)-xs(k));p(k)=(1-K(k))*ps(k);endplot(xx);holdonplot(z,'r')plot(s,'k')共128頁206共128頁207共128頁208共128頁209第三章自適應濾波系統

本章重點：掌握自適應濾波器設計的基本思想（基本原理）;掌握基於LMS演算法的自適應濾波器設計方法；瞭解熟悉自適應濾波器的應用條件及對應；掌握演算法編程的技巧。Why？維納濾波參數是固定的，適用於平穩隨機信號。卡爾曼濾波器參數是時變的，適用於非平穩隨機信號。然而，只有對信號和雜訊的統計特性先驗已知條件下，這兩種濾波器才能獲得最優濾波。What？所謂的自適應濾波，就是利用前一時刻已獲得濾波器參數等結果，自動地調節現時刻的濾波器參數，以適應信號或雜訊未知的或隨時間變化的統計特性，從而實現最優濾波。設計自適應濾波器時可以不必要求預先知道信號和雜訊的自相關函數，而且在濾波過程中信號與雜訊的自相關函數即使隨時間作慢變化它也能自動適應，自動調節到滿足最小均方差的要求（因此實際同WF及KF是一致的）。這些都是它突出優點。需要已知的內容：輸入及理想的輸出

（參考信號）Example：基於最陡下降法其中：其基本結構可表示為：應用需求（HOW？）§3.1LMS自適應濾波器基本原理介紹了自適應濾波器的結構形式、物理意義及實質，從而說明其實質是同WF一致的，只是演算法的思想不同。橫向結構的FIR濾波器的物理意義系統的零狀態回應為：語言可描述：輸出是N個所有過去各輸入的線性加權之和，其加權係數就是

上式可等價為：

即：自适应滤波器可看成是自适应线性组合器，其结构图如下：線性組合器（通用的結構模型）橫向FIR結構

橫向FIR自適應濾波器簡化形式:

自適應濾波器的工作原理介紹原則：最小均方差。需要求的是加權係數W推導的思路：將輸出y的運算式代入均方差公式中，對欲求的權係數求一階導數為零的值，即為最佳值。推導過程（驗證同KF一致性）令：則：

均方誤差的梯度（用表示）

就可得到最佳權向量，用表示，即上式是維納-霍夫方程的矩陣形式因此關鍵在於怎樣能簡便地尋找，或者說用什麼樣的演算法來求得，最常用的演算法是所謂最小均方(LeastMeanSquare)演算法，簡稱LMS演算法。其他的一些演算法(RLS遞推最小二乘)，自己看參考書。§3.2Widrow-HoffLMS演算法（遞推算法）最陡下降法原理由上面的討論得出以下結論：1：自適應濾波的結果同維納濾波器的結果一致，均滿足正交定理。2：自適應濾波器設計可採用遞推的方法求出最佳的權係數W；且可表示為演算法簡單，不涉及矩陣的複雜運算。3：該方法的關鍵技術是如何適時地求得（或估計）梯度在實際中，為了便於即時系統實現，取單個樣本誤差的平方的梯度作為均方誤差梯度的估計

上兩式的這種演算法即稱為：

Widrow-HoffLMS演算法該演算法的好處：只存在乘法和加法，簡單易於即時系統的實現。

應用Widrow-HoffLMS演算法的自適應橫向濾波器示於下圖。3.2.2能使LMS演算法收斂於的μ值範圍本節討論以下問題：1。控制參數μ是否可以任意選取，只要大於零就行嗎？2。其大小影響了演算法的那些方面？3。具體在濾波器的實際設計中，應如何確定？1。先討論LMS演算法收斂於的條件上式要收斂必須滿足：所以收斂的充分條件可以寫成2。μ是一個控制穩定性和收斂速度的參量

，其對收斂的影響有以下結論：1：在<的範圍內，當μ取得愈大，有愈大，收斂愈快。2；在滿足<<2條件下的μ，最後仍可收斂於。圖3.2.4(a)表示在這種情況收斂的過程（這種情況稱為欠阻尼情況。對應於欠阻尼情況，（1）情況成為過阻尼）。當μ選得過大，使>2時，此時將不能收斂於3.2.3LMS演算法的動態特性—學習曲線及其時間常數3.2.4梯度雜訊及其所引起的失調量[例]設M=10%（一般M=10%可以滿足大多數工程設計的要求）並設N=10，問應該取多少次迭代數？[解]按式得 所以（以迭代次數計）按經驗實際迭代次數應取100（=10

濾波器長度N）或取3~5倍。§3.3應用－－自適應抵消器§3.4自適應抵消器作為陷波器的例子

§3.5自適應濾波系統在信號處理中的其他應用3.5.1自適應仿模（Adaptivemodeling）系統3.5.2自適應逆濾波系統3.5.3參考輸入是延時k步的原始輸入的自適應抵消器

%自適應濾波程式echooff;t=(0:.01:10-0.01)';n=size(t);d=0.5*sin(2*pi*t);%參考信號noise=rand(n)-0.5;%干擾信號x=d+noise;%輸入信號M=20;%濾波器長度u=0.002;%收斂因數w=zeros(M);

fork=1:M%序列長度小於濾波器階數

y(k)=0;fori=1:k-1y(k)=y(k)+w(i)*x(k-i);ende=d(k)-y(k);fori=1:k-1w(i)=w(i)+2.0*u*e*x(k-i);endendfork=M+1:n%序列長度大於濾波器階數

y(k)=0;fori=1:My(k)=y(k)+w(i)*x(k-i);ende=d(k)-y(k);fori=1:Mw(i)=w(i)+2.0*u*e*x(k-i);endend

figure(1);plot(t,x,'b',t,y,'r',t,d,'k');xlabel('時間t');ylabel('幅值');自適應濾波實驗輸入是信噪比為1的信號，其中，v(i)是均值為零，方差為1的高斯雜訊。在這裏，我們就直接認為參考信號d(i)＝s(i)，濾波器的長度設置為20點。收斂因數u=0.005。其MATLAB根源程式如下所示：%自適應濾波的演示程式%參考的系統如下：%實際輸入信號x(n)=2^0.5*sin(0.05*pi*n)+v(n),%v(i)是均值為零，方差為1的高斯雜訊,實際輸入信號的信噪比是1%參考信號直接設為d=s;%收斂因數u=0.0005%濾波器的長度為20點，即w的長度是20點n=0:0.1:120;%系統賦初值s=2^0.5*sin(0.05*pi*n);%有用信號v=1*randn(size(n));%雜訊信號x=s+v;%實際的輸入信號d=s;%參考信號u=0.0005;%收斂因數fori=1:20,w(i)=0;%w(i)的初值都設置為0endfori=1:1200,%對y賦初值為0y(i)=0;endfori=1:1200,ifi<20forj=1:i,y(i)=y(i)+w(j)*x(i+1-j);%用卷積求y(i)ende=d(i)-y(i);forj=1:i,w(j)=w(j)+2*u*e*x(i+1-j);%修正w(n)endelseforj=1:20,y(i)=y(i)+w(j)*x(i+1-j);%用卷積求y(i)ende=d(i)-y(i);forn=1:20;w(n)=w(n)+2*u*e*x(i+1-n);%修正w(n)endendendn=1:1200;plot(n,x(n),'g',n,y(n),'r',n,s(n),'b');title('自適應濾波演示程式');xlabel('採集的點數');ylabel('各信號的幅值');legend('實際信號','輸出信號','參考信號')下麵展示了不同的u導致的濾波波形的變化u＝0.005時u＝0.01時u＝0.0005時候u＝0.05時候第五章功率譜估計

§5.1引言功率譜的物理意義及定義(維納-欽辛定理)自相關函數:對取統計平均有:計算功率譜的以上兩個公式表明:1功率譜(隨機變數)是無限個自相關函數的函數,但工程實際觀察數據只是有限個.因此分析計算隨機序列的功率譜是個功率譜的估計問題.2.不同的統計估計準則,分析估計結果是不同的.因此對功率譜估計有很多種方法.本章就是研究功率譜估計的方法.重點在於討論各種方法的特點及其實現.統計估計的一般準則通常一個好的估計，要求估計值的概率密度曲線符合正態分佈、曲線必須要窄（方差小）、且比較集中在其估計量的真值附近。評價的標準：1、偏移性B：(有偏、無偏、漸近無偏估計）定義為統計平均同真值之間的差值。2、有效性質：估計量的方差反映的是無偏估計量在真值附近的擺動情況。3、一致性：估計量的均方誤差（有效性）比較兩個有偏估計量的好壞。圖5.1二種估計的概率密度分佈５．２統計量的估計

一、均值的估計已知樣本數據：均值的估計量用下式計算（樣本之間不相關）：可以證明：當數據內部不相關時，按照上式進行估計均值是一種無偏的一致估計。方差的估計已知樣本數據：方差的估計量用下式計算（樣本之間不相關）：可以證明：這是一致估計當均值不知道時：有偏估計，但漸進無偏：無偏一致估計：隨機序列自相關函數的估計已知隨機序列的一段樣本數據，來利用這段數據估計自相關函數的方法有：１：無偏估計：該估計的特點：儘管是無偏估計，但不是一致估計，方差很大，不是一種好的估計方法２：有偏估計：該估計的特點：儘管求平均時，只用N去除不合理，但其服從漸進一致估計的原則，比無偏估計法誤差小，因此工程上常用該方法計算．利用自相關函數計算功率譜的實質圖5.2三種不同解析度的譜估計方法的例子(a)經典法BTPSD法；(b)最大熵譜估計法(自回歸PSD法)；(c)Pisavcnko諧波分解法。功率譜估計的應用

§5.3週期圖作為功率譜的估計分析該估計的性能：不難看出：w（m）是個三角函數（兩個矩形函數的卷積），被稱為Bartlett窗函數，用表示。自相關函數的真值當時說明：１：方差不等於零，不滿足一致估計的條件．２：方差較大，不是功率譜的好估計．解決辦法：將週期圖進行平滑（或平均）處理，滿足一致性．圖5.3高斯白雜訊序列的週期圖

§5.4平滑後的週期圖作為PSD的估計

5.4.1巴特利特(Bartlett)平均週期圖的方法[例如]為了說明經平均後的週期圖作為功率譜估計的實際效果，設有一零均值高斯分佈的隨機過程，其功率譜密度為圖5.4與的特性

圖5.5平滑後的週期圖（每段取8個數據）

圖5.6平均後的週期圖(每段取16個數據)

5.4.2窗口處理法平滑週期圖

5.4.3Welch法主要針對Bartlett法提出二方面的修正：其一是選擇適當的窗函數，並在週期圖計算前直接加進，這樣保證了對於任何的窗函數，功率譜密度估計非負。其二是分段之間可以有重疊（相當於數據長度變長了分段多了）因此方差減小。經典譜估計法總結：實質：根據測得的數據，通過計算FFT（傅立葉變換）獲得相應的譜密度資訊。特點：計算簡單，效率高；缺點：由於將測得數據以外的數據均看作零（實質相當於乘上了一個矩形窗口），因此頻率解析度低。現代譜估計法的基本思想：克服經典譜估計的缺點，另辟新徑，從數學模型（參數）的角度出發。認為已觀察到的數據是白雜訊通過一個數學模型產生的。處理思想及步驟：1選擇一個合適的模型

2確定模型參數（根據觀察值）

3由模型求PSD

關鍵1：模型選擇問題（AR，MA，ARMA）

2：參數確定方法（導致產生了各種演算法）

§5.5自回歸模型法(AR)數字系統的數學模型:差分方程求功率譜的實質變為確定系統參數的問題ARMA模型圖5.7自回歸模型到底選擇什麼模型？為何討論AR模型？1。Wold分解定理：任何有限方差的ARMA或MA平移過程可以用可能是無限階的AR模型表達；同樣，任何ARMA或AR模型可以用可能是無限階的MA模型表達。因此，如果在這三個模型中選了一個與信號不匹配的模型，利用高的階數仍然可以得到好的逼近。2。由於對AR模型參數的估計，得到的是線性方程。故AR模型比ARMA以及MA模型有在計算上的優點。同時，實際的物理系統往往是全極點系統。所以研究有理分式傳遞函數的模型，主要研究AR模型。如何根據自相關函數確定系統參數最小均方誤差準則下的線性預測器求AR模型各參數的另一種方法：自適應濾波器原理實質上：建模和預測是等價的。

圖5.8AR模型的H(z)與其逆濾波器的串接圖5.9橫向結構的預測誤差濾波器

§5.6最大熵譜估計法（MESE）

MaximumEntropySpectralEstimation）

5.6.1最大熵譜估計法的基本思想及其與線性預測AR模型法的等價性熵的基本概念：熵代表一種不定度，最大熵代表最隨機，意味著對應的PSD最平坦。其定義為：結論：

上式實質為Yuler－walker方程，因此其同AR模型等價

5.6.2Levinson-Durbin遞推算法：

Yule－Walker方程的一種高效解

§5.7預測誤差格型濾波器及伯格(Burg)遞推算法

5.7.1預測誤差格型濾波器圖5.10格型預測誤差濾波器圖5.11後向預測誤差

5.7.2Burg遞推算法——Kp的確定

5.7.3正反向線性預測最小二乘法第九章 小波變換及其應用

內容提要：本章主要由傅裏葉變換的不足引出小波變換的概念和特點，介紹小波變換的基本理論和分析步驟，並以MATLAB為工具介紹小波分析的基本命令和應用實例。第一節從傅裏葉變換到小波變換前面重點講述了確定性信號的頻域分析、應用以及隨機信號分析的基本內容。本節開始介紹常用於非穩態分析的一種新的變換方法——小波變換。傅裏葉變換和傅裏葉反變換是同一能量信號的兩種不同表現形式，正變換把一個信號函數分解成眾多的頻率成份，反變換將這些頻率成份重構為原來的信號，轉換過程中能量保持不變。在採樣定理的指導下，離散傅裏葉變換（DFT）為電腦實現傅裏葉變換提供了理論基礎，快速傅裏葉變換（FFT）的提出則使傅裏葉分析變得實用。然而，傅裏葉變換的作用是有限的，這裏主要談兩點。一是傅裏葉變換沒有包含時間資訊。實際信號中含有大量非穩態的成份，如突變、偏移、趨勢、事件的開始、終止等等，其頻率特性將隨時間而發生變化，對這些信號進行分析時需要提取某一時間段的頻域資訊或某一頻率段的時間資訊。二是傅裏葉變換不利於分析時變信號。時變信號可以分段進行研究，對變化快的信號，其頻率高，取小的時間間隔有利於提高分析的精度；對變化慢的低頻信號，取大的時間間隔才可以收集到完整的資訊以進行分析。而傅裏葉變換不能實現這樣的時－頻局部化分析。(a)Haar小波 (b)Morlet小波小波的波形可以進行拉伸和壓縮，以用以分析信號的輪廓和信號的細節，這種拉伸和壓縮變換被稱為尺度變換。圖9-2為不同尺度的Daubechies小波的示意圖，從圖中可以看出，尺度值a越小，圖形壓縮越厲害；反之，尺度值越大，圖形拉伸越厲害。圖9-2不同尺度小波示例圖9-3為短時傅裏葉變換和小波變換的對比示意圖。圖（a）表示短時傅立葉變換，對任意頻率，其窗口都一樣；圖（b）表示小波變換，它將一個時間信號轉換為時間和頻率的二維函數，可以提供信號在某個時間段和某個頻率範圍的一定資訊。STFT圖9-3短時傅立葉變換和小波變換示意圖WFT第二節小波和小波變換一、母小波和小波1．定義滿足如下允許條件或其等價條件的函數稱為一個母小波函數（MotherWaveletFunction），其中為的傅裏葉變換。式（9－2）說明母小波函數具有一定的振盪性，即包含某種頻率特性。對母小波函數作伸縮、平移得式中，，稱為小波函數，簡稱小波。式中變數a反映函數的尺度（或寬度），a越小，則波形壓縮越厲害，a越大，則波形拉伸越多，變數b檢測沿t軸的平移位置。一般情況下，母小波函數能量集中在原點，小波函數能量集中在b點。母小波具有以下三條性質：①，即母小波具有零直流分量。②母小波及其生成的小波函數均為帶通信號。③母小波及其生成的小波函數隨t的延伸而快速衰減。2．小波的尺度、時－頻關係及濾波特性設母小波函數的窗口寬度為，窗口中心為，則相應的連續小波的窗口中心為窗口寬度為即時、頻窗口的面積與尺度a無關，它的與的大小是相互制約的，這正是海森堡測不准原理的要旨，即時間解析度和頻率解析度是相互制約的。由此，可以得到以下結論：①小波的尺度a和頻率相對應。小尺度a對應高頻的，大尺度a對應低頻的。②時間、頻率（尺度）窗口的形狀隨a而變化。這是和短時傅裏葉變換的等時間窗所不同的。③不可能同時得到很高的時間、尺度解析度。這是由於時間、尺度窗口的面積恒定決定的。④小波可以看作是一組衡Q濾波器。由於小波母函數具有有限頻帶，其伸縮、平移得到的小波為一組帶通濾波器。帶通濾波器的品質因素被定義為通帶寬度與中心頻率的比值，設母函數的品質因素則由它生成的小波的品質因素不隨尺度a而變化。二、連續小波變換和反變換1．連續小波變換連續小波變換有兩個參數a和b，對同一個a，信號可以分解成不同時移b的小波的疊加，而改變a值，同一個信號又可以在不同層次上進行由粗到細的分解，獲得小波變換的步驟可以概括為以下五步：選擇尺度a一定的小波，將它與被分析信號的初始段進行比較，如圖9-4（a）。通過計算得出該段信號與所選小波的相關程度，越大，說明二者越相似。這一結果依賴於所選擇小波的形狀。改變b值實現小波平移，如圖9-4（b），重複（1）、（2）步，直至對信號完成一次比較分析。增加尺度因數a，即拉伸小波，重複（1）—（3）步，對信號進行下一輪分析，如圖9-4（c）。對全部尺度因數重複（1）－（4）步，可以得到使用不同尺度評估信號在不同時間段的大量係數值。這些係數反映了被分析信號在小波函數上的投影，可以用灰度圖表示出來。圖9-4（d）所示是含有隨機雜訊的正弦波經過小波變化後得到的灰度圖，從圖中可以看出，尺度大時可以提取出信號的輪廓，將高頻的雜訊去除，從而發現被分析信號具有明顯的週期性。（a） （b） 圖9-4小波分析的步驟2．小波變換的性質①線性②時移不變性設的小波變換為，則的小波變換為。③尺度變換特性設的小波變換為，則的小波變換為。④微分特性⑤能量守恆特性⑥冗餘度連續小波變換把一維信號變換到二維空間，因此小波變換中存在多餘的資訊，即為冗餘度。度量冗餘度的量稱為再生核，它反映了小波變換二維空間兩點之間的相關性。3．連續小波反變換此式說明了根據精確恢復信號的方法。三、離散小波變換離散小波變換中，取小的參數a、b將有助於提高重構信號的精度。四、二進小波變換和正交小波變換若連續小波僅對尺度量進行離散化，且取，則為二進小波，相應的小波變換為二進小波變換。它的特點在於尺度參量離散，而時域上的平移量仍保持連續變化。因此，與離散小波變換不同，二進小波變換保持了連續小波變換所具有的時移不變性，這是二進小波最大的特點。第三節多解析度分析與小波包分析一、多解析度分析L2(R)空間的多解析度分析是指滿足下列條件的一個空間序列：單調性：；漸進完全性：；伸縮性：對任意，有平移不變性：若，則Riesz基存在性：存在函數，使構成的Riesz基，即對任意的，存在唯一的序列，使得。V0W0V1W1W2V2L2(R)圖9-5尺度空間與小波空間的關係總結以上分析，多解析度分析的過程可以簡述如下：由條件（5）中的Reisz基可以生成尺度函數，由此可以構造出空間的標準正交基。在此基礎上，可以得到與空間正交互補的小波空間，並可建立小波空間的標準正交基。實際信號所處空間的標準正交基可以在及的標準正交基的基礎上得到，可以證明，它與二進正交小波是一致的。因此，可以由尺度函數構造出小波函數，不過設計滿意的尺度函數是比較複雜的。二、兩尺度方程

因此，為低通濾波器，為高通濾波器，分別對應於尺度函數的低通特性和小波函數的帶通特性。圖9-6多解析度分析圖示三、小波包圖9-7小波包圖示第四節小波變換的快速演算法本節介紹基本的一維Mallat演算法，Mallat演算法是一種基2類小波快速演算法，在語音和圖像處理中應用很多。此即為小波變換係數的重建公式

圖9-8Mallat演算法分解與重構演算法圖（a）分解快速演算法示意圖（b）重構快速演算法示意圖圖9-9二次分解和二次重構示意圖

（a）二次分解（b）二次重構第五節以MATLAB實現小波分析及應用一、常用小波

1．Haar小波2．Daubechies小波

Daubechies小波一般簡寫為dbN，N為小波的階數。當N=1時，為Haar小波，時，dbN沒有顯式運算式，圖9-2分別顯示了db3、db5、db7、db10的情況。3．Mexicohat小波例9-1．作出[-5,5]上取1000點的Mexicohat小波。解％作出Mexicohat小波的圖形lb=-5;ub=5;n=1000;[psi,x]=mexihat(lb,ub,n);plot(x,psi);title(‘Mexicanhatwavelet’);4．Morlet小波此小波的圖形如圖9-1（b），定義為由於複值Morlet能提取信號的幅值和相位資訊，較常用。5．Meyer小波圖9-11Meyer小波二、常用信號的小波分析1．命令行方式小波分析

圖9-12加雜訊正弦信號的小波分析

例9-2．用命令行方式對加入白雜訊的正弦信號進行連續小波分析。結果如圖9-12所示。解t=1:1000;y=sin(0.1*t);ynois=sin(0.1*t)+0.5*randn(size(t));c=cwt(y2,1:64,'db4','plot');2．圖形化方式小波分析圖9-13圖形化方式小波分析的主菜單

例9-3．以圖形化方式對例9-2進行小波分析，要求以load調用信號。解％先建立待分析信號，並以save命令保存。fort=1:1000y(t)=sin(0.1*t)+0.5*randn(size(t));endsavey;%===========================％打開小波分析工具箱主菜單wavemenu;%===========================點擊ContinuousWavelet1-D，並選擇File=>LoadSignal=>y，可調用含白雜訊的正弦信號，在生成的介面中選擇Wavelet：db4；Scale：1:1:64，點擊Analyse鍵，可得分析結果如圖9-14。圖9-14圖形化連續小波分析(a)尺度為10時的小波係數圖

圖9-15不同尺度對應的小波係數圖(b)尺度為50時小波係數圖

3．信號突變點檢測小波變換可以容易地進行信號的特徵提取。例9-4．以MATLAB自身帶的nearbrk信號進行信號突變點檢測。解%裝入原始信號loadnearbrk;s=nearbrk;%===========================%畫出該信號的波形圖subplot(4,2,1);plot(s);ylabel('s');title('原始信號s和信號的近似a、細節d');%===========================%用小波db2進行3層分解[c,l]=wavedec(s,3,'db2');fori=1:3decmpa=wrcoef('a',c,l,'db2',4-i);subplot(4,2,2*i+1);plot(decmpa);Ylabel(['a',num2str(4-i)]);endfori=1:3decmpd=wrcoef('d',c,l,'db2',4-i);subplot(4,2,2*(i+1));plot(decmpd);Ylabel(['d',num2str(4-i)]);end圖9-16例9-4運行結果4．多解析度分解與重構例9-5．給定正弦信號，對其進行多解析度分解與重構。解％生成原始信號t=0:pi/100:4*pi;s=sin(t+3*pi/4);subplot(3,2,1);plot(s);ylabel('s');gtext('原始信號');%====================================%對s進行小波分解：db13層[c,l]=wavedec(s,3,'db1');%====================================%提取小波分解的低頻係數a3a3=appcoef(c,l,'db1',3);%====================================%提取小波分解的各層高頻係數d3=detcoef(c,l,3);d2=detcoef(c,l,2);d1=detcoef(c,l,1); %====================================%繪出各係數的圖形subplot(3,2,3);plot(a3);ylabel('a3');subplot(3,2,2);plot(d3);ylabel('d3');subplot(3,2,4);plot(d2);ylabel('d2');subplot(3,2,6);plot(d1);ylabel('d1');%====================================%重構信號ss1=waverec(c,l,'db1');subplot(3,2,5);plot(s1);ylabel('s1');gtext('重構信號');圖9-17例9-5運行結果三、信號處理1．信號的自相似性檢測例9-6．利用小波分析來檢測MATLAB自帶的信號vonkoch。解loadvonkoch;s=vonkoch;subplot(2,1,1);plot(s);title('原始信號');subplot(2,1,2);f=cwt(s,[2:2:128],'coif3','plot');title('小波分解自相似指數圖');xlabel('時間');ylabel('變換尺度');圖9-18例9-6運行結果2．識別信號中的頻率成分例9-7．以小波分析識別MTALAB自帶的sumsin信號，此信號由三種不同頻率正弦信號疊加而成。解loadsumsin;s=s