BIOINGEGNERIA

S. Salinari

Lezione 2

Filtro di Wiener Stimatore ottimo non recorsivo

Il filtro di Wiener è lo stimatore ottimo non recorsivo dal punto di vista dell’errore quadratico medio pe = E(e2).

IPOTESI: 1. y(k) = x + v(k) non essenziale

2. x e v processi stazionari

3. M

Il filtro non lavora il tempo reale. Gli m campioni di y sono tutti acquisiti

Si ha un filtro lineare tempo-invariante per cui y* x* =i h(i) y(i) i=1...M

le h(i) vanno calcolate minimizzando pe = E(e2)

pe = E[x - i h(i) y(i) ]2

pe/ h(j) = -2 E[x - i h(i) y(i) ]y(j) = 0 j=1...M E[e y(j)] = 0

i h(i) E[y(i) y(j)] = E[xy(j)] j=1...M i h(i) py(i,j) = pxy(j)

py(i,j) e pxy(j) noti, h(i) incognita

pe =E(e2)=E{e[ x - i h(i) y(i) ]}=[ E(x2) - i h(i) E[xy(i)] = E(x2) - i h(i) pxy(i)

Tenendo presente che la matrice Py è simmetrica py(i,j) = py(j,i)

In forma matriciale si ottiene: h = Py-1 pxy EQUAZIONE DI WIENER-HOPF

x*= hT y =pxyT Py

-1 y

pe= E(x2) – pxyT Py

-1 pxy

Con h, y vettori colonna (M1)

Py matrice di autocorrelazione (M M)

pxy vettore di mutua correlazione (M1)

Osservazioni pratiche: Se ci si discosta dalle ipotesi di base si ha un filtro subottimo.

Filtro di Wiener Stimatore ottimo non recorsivo

Filtro di Wiener Stimatore ottimo non recorsivo

APPLICAZIONI – Esempio 1

v2 se k=j

y(k) = x + v(k) con: E(x) =0; E(x2) = x2; E(v)=0; E[v(k)v(j)] =

0 se kj

La soluzione con la formula di Wiener-Hopf porta a:

h(1) = h(2) = ...........= h(M) =

x* = =

pe = per M grande pe

22

2

vx

x

M

M

)i(yM

i 1

2

2

x

v

Mv2

Mv 2

Filtro di Wiener Stimatore ottimo non recorsivo

APPLICAZIONI – Esempio 2

y(k) = x ·k + v(k) rampa con coefficiente angolare x

E(x) =x0; E(x2) = [E(x)]2+ x2 =S; E(v)=0; E[v(k)v(j)] = v

2 per i=k

0 per ik

Stimare la x* noti 2 punti per k=1 e k=2

(S+v2)h(1) + 2Sh(2) = S

2Sh(1) + (4S+v2)h(2) = 2S

jS)x(jE)]j(vx(x[E)]j(xy[E)j(p

)j,i(ijS)j,i()x(ijE

)]}j(vjx)][i(vix{[E)]j(y)i(y[E)j,i(p

xy

vv

y

2

222

SSp

S

)(y)(y*x

)(h)(h

e

x

v

221

5

22

5

11

2

2

Filtro di Wiener Stimatore ottimo non recorsivo

SVANTAGGI DELLA FORMULAZIONE NON RECORSIVA

I principali svantaggi nell’applicazione della formulazione non recorsiva del filtro di Wiener sono:

1. Richiede la conoscenza a priori della Py e della pxy.

2. Il numero di campioni m viene definito a priori

3. Se m cambia per qualunque ragione (ad es. più dati disponibili) il calcolo va ripetuto

4. Richiede l’inversione di una matrice (mxm) Py che può essere laboriosa.

Filtro di Wiener Formulazione recorsiva

)j,i()]j(v)k(v[E;)v(E;)x(E;)x(E)k(vx)k(y vx 222 00

Nel caso in cui:

Si può ottenere una formulazione recorsiva ricordando che: h(k) =

Infatti si può scrivere :

kk vx

x 122

2

)]k(y)k(*x)[k(a)k(y)k(*x

-a(k)b(k)

)k(y)k(b)k(*x)k(a

)k(yk

)k(*xk

k

k

)k(k

k

)i(y

k

k)x*(k

)k(y)i(yk

)k(*x

k

i

k

i

111

:ha si 1 Poichè

1

11

1

11

1

11

ottiene sik per dividendo e ndomoltiplica

11

11

1

1

)j,i())j(v)i(v(E;)v(E;)x(E;)x(E)k(vx)k(y vx 222 0 0con

Ta(k)

y(k+1) x*(k+1)

x*(k)

++

+-

Filtro di Wiener Segnale variabile nel tempo

Consideriamo che y(k) = x(k)+v(k) con y(k) e x(k) campioni di grandezze scalari x(t) e y(t).

La misura varia quindi non solo per il contributo del rumore ma anche per la variabilità nel tempo del segnale.

Minimizzando l’errore quadratico medio per ogni k si otterrà un’espressione analoga al caso in cui x era costante e precisamente:

pe(k) = E[x(k) - i h(i) y(i) ]2

pe(k) / h(j) = -2 E[x(k) - i h(i) y(i) ]y(j) = 0 j=1...M

i h(i) E[y(i) y(j)] = E[x(k)y(j)] j=1...M

Poichè si dovrà scrivere un’equazione del tipo per ogni k si avrà:

i h(k,i) E[y(i) y(j)] = E[x(k)y(j)] j,k=1...M

In forma vettoriale si avrà:

HE(yyT) = E(xyT)

h(1,1)...h(1,j)...h(1,m) E[x(1)y(1)]...E(x(1)y(m) E[y(1)y(1)]...E(y(1)y(m)

con H = h(k,1)...h(k,j)...h(k,m E(xyT) = E[x(k)y(1)]...E(x(k)y(m) E(y yT) = E[y(k)y(1)]...E(y(k)y(m)

h(m,1)...h(m,j)...h(m,m) E[x(m)y(1)]...E(x(m)y(m) E[y(m)y(1)]...E(y(m)y(m)

Filtro di Wiener Segnale e Osservazione vettoriale

Si consideri ora il caso in cui, in k, si abbiano q segnali x: x1(k), x2(k),...xq(k) ed r osservazioni y: y1(k), y2(k),...yr(k).

In tal caso si ha: H E(YYT)=E(xYT) H Py=Pxy

H(1,1) H(1,2).....H(1,m) hki(1,1)...hki(1,r)

H = ................................... H(k,i) = .............

H(m,1) H(m,2).....H(m,m) hki(q,1)...hki(q,r)

Py(1,1) Py(1,2)..... Py(1,m) pyij(1,1)... pyij(1,r)

Py = ................................... Py(i,j) = .............

Py(m,1) Py(m,2)..... Py(m,m) pyij(r,1)... pyij(r,r)

Pxy(1,1) Pxy(1,2)..... Pxy(1,m) pxyij(1,1)... pxyij(1,r)

Pxy = ................................... Pxy(i,j) = .............

Pxy(m,1) Pxy(m,2)..... Pxy(m,m) pxyij(q,1)... pxyij(q,r)