第二十八讲参数化功率谱估计方法

第二十八讲《现代数字信号处理》江

汉：理工大学国家短波工程技术理工大学通信第14章：功率谱估计第二十六讲：14.1、(

著)3.6

节基于有限长信号观察的功率谱估计（掌握）（周期图法、离散

变换DFT）第二十七讲：14.2

节功率谱估计的非参数化方法（掌握）（改进的周期图法）第二十八讲：14.3

～

14.4节功率谱估计的参数化方法（掌握）（AR、MA、ARMA模型）滤波器组方法（掌握）第二十九讲：课堂 ：功率谱估计的仿真理工大学通信回顾：非参数化方法

-

周期图理工大学通信xx1KK

1P

(

f

)

i021x

(n)eMM

1

j

2π

fnn0

ii

P

(

f

)iP

(f)

xxxxWelch方法：21xxMUx

(n)w(n)eM

1

j

2π

fnn

0

iPi

(

f

)

ixxP

(

f

)1LL

1i

0Pxx

(

f

)

Q

8L

16N9

9MwMQ

K

N

B质量因子：

1QPBartlett方法：

分段+平均平稳随机信号:x(n)

n

0,1,...,

N

1

N

K

M2xx1NN

1n

0周期图法：P

(

f

)

x(n)e

j

2π

fnxi

(n)

x(n

i

M

)+加窗ix

(n)

x(n

i

D)Blackman&Turkey方法：BTQ

3N2M自相关加窗M

1

j

2π

fmPxx

(

f

)

rxx

(m)

w(m)en

M

1回顾：非参数化方法

-

周期图理工大学通信四种方法回顾：非参数化方法

-

周期图特点非参数化方法的特点：基于序列或自相关的DFT,计算简单；高分辨率需要较长数据，加窗存在频谱泄漏；功率谱估计性能（质量因子）逐步提高。周期图法的共性分析：

假设合理性？ 能否去掉？从DFT角度看，功率谱密度估计是周期的，隐含假设数据也是周期的（周期为N）；从自相关函数角度看，只有2N-1点，隐含假设当|m|≥N时，自相关值为0。理工大学通信回顾：非参数化方法

-

周期图14.3

功率谱估计的参数化方法理工大学通信参数化基本思想114.3

功率谱估计的参数化方法理工大学通信白噪声参数化模型：将数据序列x(n)建模为一个线性系统的输出。有理系统函数：参数化基本思想1参数化？qpb

zB(z)A(z)

kk

0

kk

1H

(z)

k1

ak

zp

qx(n)

ak

x(n

k

)

bk

w(n

k

)k

1

k

0w(n)x(n)H

(z)差分方程：平稳随机过程理工大学通信1）输入信号为白噪声序列，功率谱密度为 ，通过滤2wσ波器

H(f)后，输出为平稳随机过程

xx(

f

)

；2）输入信号为平稳随机过程，功率谱密度为

xx

(f)，通过滤波器1/H(f)后，输出为白噪声序列。w(n)H

(z)w(n)1

/

H

(z)2 (

f

)

σ

2

H

(

f

)xx

w参数化基本思想1平稳随机过程wB(

f

)

2A(

f

)

2

σ

2

功率谱密度：两步：先估计模型参数、后计算功率谱密度；x(n)Wold表示x(n)理工大学通信引入自相关?*xxγ

(m)

Ex(n)

x

(n

m)参数化基本思想1自相关与模型参数xxγ

(m)

k

xxa

γ

(m

k

)

bkpk

1p

qx(n)

ak

x(n

k

)

bk

w(n

k

)k

1

k

0左侧只有x(n):

ak=0,k>=1左、右侧均有延迟项：1) 右侧只有w(n):

b0=1,

bk=0,k>0

自回归（AR）滑动平均（MA）ARMA求

bk

,

ak

?2wq

σ

h(k

m)k

0理工大学通信ppxxk

1q

m2k

0k

1m

q0

m

qxxγ

(m

)

γ

*

(

m)m

0

ak

γxx

(m

k

)k

xxa

γ

(m

k

)

σw

h(k

)bk

mARMA(p,q)过程：参数化基本思想1ARMA模型参数xxxxxxxxxxxxγ

(q

1)γ

(q)γxx

(q

p

1)

a1

γxx

(q

1)

γγ

(q

p+1)

a

(q

1)

2

γ

a

γ(q

p

1)γ

(q

p

2)(q

p)

xx

xx

p

γ

(q)γxx

(q)γxx

(q-1)

线性模型扩展了自相关函数m范围。Γ

p

a

p

γ

pppp矩阵形式：

a

Γ1

γ

理工大学通信xxxxγ

(0)γxx

(

p)

1

σ

2

w

γxx

(0)

γxx

(1)γxx

(

p+1)

a1

0

γ

a

0

(

p)

γ

(

p

1)

xx

p

γxx

(-1)

γxx

(0)

AR（p）过程：pa

γ

(

m

k

)

σ

2

m

0k

xx

wγxx

*

(

m

)m

0xxγ

(

m

)

m

0

ak

γxx

(

m

k

)k

1pk

1AR模型参数Yule-Walker(尤尔-

)方程Toeplitz矩阵形式：

系数ak+

rxx(m≤p)

->rxx(m>p)参数化基本思想1理工大学通信参数化基本思想1AR模型参数求解214.3

功率谱估计的参数化方法理工大学通信xxxxσ

2

γxx

(

p)

1

w

γxx

(0)

γxx

(1)γxx

(

p+1)

a1

0

γ

a

0

xx

p

(

p)

γ

(

p

1)

γxx

(-1)

γxx

(0)

γ

(0)

γxx

(1)

γxx

(0)

γxx

(

p

1)γxx

(

p)ap

pppΓ

a

γL-D方法p

γxx

(

p

1)

a1

γxx

(0)pa

Γ1

γpY-W方程求解AR模型参数求解2AR（p）的Yule-Walker方程矩阵形式：理工大学通信pa

Γ1

γp

pAR模型参数估计2自相关计算：*xxNN

m1n0r

(m)

1x(n

m)x

(n)*xxxxr

(m)

r

(m)

有限长:*xxγ

(m)

E

x(n

m)x

(n)m

0,1,

2,

N

1m

0,1,

2,

p

12xxw(

f

)

σ

2

H

(

f

)w

σ

2

B(

f

)

2A(

f

)

21wA(

f

)

2

σ

2

功率谱密度：

[ap

(0),

ap

(1)

,...,

ap

(

p

1)]2wpxxppaˆ

(k

)eσˆ

2

j

2π

fkPYW

(

f

)1

k

1pxxwp

xx

pk

1(k

)σˆ

2

γ

(0)

a

(k

)γMSE:Y-W方程求解理工大学通信阶数p的选择？AR模型参数求解22wpxxpaˆ

(k

)eσˆ

2

j

2π

fkk

1PYW

(

f

)1

p2pwpk

1σ

γxx

(0)

ap

(k

)γxx

(k

)阶数p的常用准则：Y-W方程求解：计算顺序：系数、MSE、功率谱密度Akaike信息准则（AIC）：最小描述长度（MDL）

：

2

p

/

NwpMDL(

p)

N

ln

σˆ

2

p

/N修正，准确性优于AIC；在某个p值获得最小值；数据长度N较小时，阶数：N/3~

N/2理工大学通信其它求解？AR模型参数求解2Y-W方程求解：最小均方误差准则（MMSE）fpf

p

(n)2

ε

E等价于前向 误差2n

mN

12

εm

mf

(n)

gm

(n)Burg求解：

基于线性的最小二乘方准则（LS）约束条件：AR参数满足L-D递归非约束最小二乘（LS）求解：2xxpEˆ

paˆ

(k

)e

j

2π

fkPBU

(

f

)

1

pk

1xx

2paˆ

(k

)eε

LS

j

2π

fkk

1PLS

(f

)

p

1

p2wpxxpaˆ

(k

)eσˆ

2

j

2π

fkk

1PYW

(

f

)

1

p理工大学通信AR与ARMA模型关系AR模型参数求解2AR（p）过程：

平稳随机过程的功率谱估计B(

f

)

2 (

f

)

σ

2

H

(

f

)

2

σ

2

xx

w

wA(

f

)

21wA(

f

)

2

σ

2

AR(p)平稳随机过程 +

噪声：1A(

f

)

2

σ

2nxx (

f

)

σ

2

wxx2nA(

z)

2σ

2

(

z)

w

σnσ

2A(z)A(z

1

)

w

σ

2σ

2

σ

2

A(

z)

A(

z

1

)

w

n

A(

z)

A(z

1

)等价于ARMA(p,p)过程在噪声环境下，ARMA(p,p)模型参数估计更好。理工大学通信性能比较AR模型参数求解2△f=0.13AR（4）模型理工大学通信AR模型参数求解2△f=0.07性能比较AR（4）模型理工大学通信AR模型参数求解2性能比较ARMA（10，10）模型适合于噪声环境理工大学通信回顾：非参数化方法

-

周期图14.3

功率谱估计的参数化方法14.4

功率谱估计的滤波器组方法理工大学通信功率谱密度的滤波器组表示：14.4

功率谱估计的滤波器组方法特定频率的功率谱估计整个功率谱密度的估计由一组FIR滤波器来实现。如何设计滤波器组？理工大学通信周期图的滤波器组实现114.4

功率谱估计的滤波器组方法理工大学通信周期图法：x(n)

n

0,1,...,

N

1平稳随机信号:21xxNN

1n0P

(

f

)

x(n)e

j

2π

fn周期图法21xx

kNN

1n0P

(

f

)

x(n)ek

j

2πf

nk给定频率:

fN

f21NN

1n

0x(n)ekj

2π

f

(

N

n

)定义FIR滤波器:1kNe

j

2π

fk

nh

(n)

n

02N

1n

0Pxx(

fk

)

x(n)hk

(

N

n)2N

1n

0l

N

x(n)hk

(l

n)1周期图的滤波器组实现理工大学通信对应FIR滤波器:1e

j

2π

fk

nkNh

(n)

n

0

fNfkkkh

(n)z

nn

0H

(z)

1/

N

1

e

j

2π

fk

z1单极点（单位圆上）2N

1n

0Pxx(

fk

)

x(n)hk

(

N

n)滤波器组周期图的滤波器组实现1滤波器组理工大学通信1对应FIR滤波器:hk

(n)Nn

0

~

N

1fkN

f矩形窗21xxNN

1P

(

f

)

n0x(n)e

j

2π

fn21NN

1n0

x(n)e

j

2π

fn

2xxkx(n)e1NN

1

j

2π

fn

h

(

N

n)P

(

f

)

n0其它窗?更好的滤波器组?周期图的滤波器组实现1其它滤波器组？周期图的滤波器组固定，没有考虑数据统计特征。理工大学通信周期图的滤波器组实现1最小方差（MV）谱估计214.4

功率谱估计的滤波器组方法理工大学通信x(n)

n

0,1,...,

N

1平稳随机信号:y(n)

h(k

)x(n

k

)设计FIR滤波器：

h(k

)

k

0

~

pp最小方差（MV）谱估计2输出响应：k

0Xt

(n)

[

x(n)

x(n)

Xt

(n)

h

x(n

p)]输出方差：σ22yH

*t

Ey(n)

E

h X

(n)

X

(n)h

h

H

Γ

hxx通过FIR滤波器后输出:自相关矩阵pk

0

1h(k

)e

j

2πkfl频率e

j

2π

fl约束条件：lEH

(

f

)

h

1...le

j

2π

fle

j

2π

pfl

EH

(

f

)

1基本思想理工大学通信最小方差滤波器:目标：设计的滤波器能够无失真传递频率分量，而远离频率

的分量严重衰减。最小方差（MV）谱估计22H

tHxx2yy(n)σ

E

E

h X(n)

X

(n)h

h

Γhfl约束条件：EH

(

f

)

h

1lMV准则fl模型求解：

日乘子法J

hH

Γh

Re

λ

(E*

H(

f

)

h

1)xxl求导：xx

lJ

2Γ h

λ*

E(

f

)

0*2xx

l1hh

λ

Γ

E(

f

)

λ*

2EH

(