Filtro di Wiener

Il Filtro di Wiener è un filtro per l'elaborazione dei segnali su base statistica, proposto da Norbert Wiener negli anni 1940 e pubblicato nel 1949.

Il filtro viene oggi realizzato con tecniche digitali, ma può anche essere realizzato con sistemi analogici, come il prototipo costruito al MIT sulla base del progetto di Wiener.

In generale, l'obiettivo di un filtro è l'eliminazione del rumore presente in un segnale.

Tipicamente i filtri sono concepiti per una ben precisa risposta in frequenza e possono operare la separazione tra segnale e rumore a condizione che questi occupino diverse bande di frequenza. Il filtro di Wiener supera questa limitazione affrontando il problema del filtraggio con un approccio statistico. Si assume di avere conoscenza delle caratteristiche spettrali del segnale originale e del rumore, e si cerca il filtro LTI la cui uscita sia "il più vicina possibile" al segnale originale, indicando con questa espressione la minimizzazione di una prefissata misura dell'errore compiuto nell'operazione.

I filtri di Wiener sono caratterizzati dalle seguenti caratteristiche Ferreriane:

1. Assunzioni: il segnale e il rumore additivo, ovvero che si somma al segnale, sono processi stocastici stazionari e lineari con caratteristiche spettrali note, o con autocorrelazione e mutua correlazione note.
2. Condizioni: il filtro deve essere fisicamente realizzabile, ovvero la sua risposta (effetto) non può precedere temporalmente l'ingresso (causa). Essendo l'effetto successivo alla causa questi filtri vengono anche chiamati causali. La condizione di causalità può essere abbandonata; in questo caso si parla di filtri, o di soluzioni, non causali.
3. Criteri prestazionali: minimizzazione dell'errore quadratico medio.

Questo filtro è usato frequentemente nel processo di deconvoluzione; per questo utilizzo vedi deconvoluzione di Wiener.

Il filtro venne proposto da Norbert Wiener durante gli anni della seconda guerra mondiale, per rispondere alla necessità del governo USA di migliorare le capacità di tiro di un sistema di puntamento antiareo, basato sull'uso del radar. La velocità dei bersagli poneva problemi nuovi: il sistema doveva necessariamente essere automatizzato, e doveva essere in grado di dirigere il tiro verso la posizione futura dell'aereo al momento dell'impatto, che doveva quindi essere stimata. Era perciò necessario integrare due diverse tecnologie, quella dei controlli automatici e quella delle telecomunicazioni, fino ad allora applicate in ambiti diversi.

L'analisi dei sistemi di puntamento automatico portarono Wiener a constatare la somiglianza, in termini funzionali, del loro comportamento con quelli tipici degli esseri viventi, nello svolgimento di attività apparentemente rivolte al raggiungimento di uno scopo. Questi stimoli suggerirono a Wiener la definizione di una nuova scienza, la cibernetica, che studiasse da un punto di vista matematico unitario comportamenti comuni a sistemi sia naturali che artificiali^[1].

Il sistema di puntamento doveva elaborare del segnale radar ricevuto, al fine di eliminare il rumore in esso contenuto, e stimare la posizione futura del bersaglio. Alla classica operazione di filtraggio del segnale (che non poteva però essere basata, in questo caso, sulla separazione delle frequenze) si aggiungeva quindi una nuova esigenza di predizione del suo futuro andamento, sulla base della sua storia passata. Questi due casi, apparentemente distinti, vennerono compresi da Wiener all'interno di un'unica teoria, che egli sviluppò riconducendo il problema ad una classe di equazioni integrali da lui già studiate nel 1931, in collaborazione con il matematico ed astrofisico austriaco Eberhardt Hopf (1902 - 1983)^[2][3]. Il primo documento che descrive il filtro è del 1942, ma poiché Wiener lo aveva sviluppato per una applicazione militare, esso fu secretato e non poté circolare fino al 1949, quando ne fu permessa la divulgazione^[1].

Nel 1941 Andrey Kolmogorov aveva pubblicato un lavoro equivalente sul problema della predizione, nel caso a tempo discreto^[4], successivamente esteso al caso tempo continuo dal matematico russo M. Krejn a partire dal 1944. A sua volta, Kolmogorov si era basato su alcuni risultati dello statistico svedese H. Wold^[5]. Tuttavia, Wiener aveva lavorato indipendentemente da questi precedenti, che non conosceva. Per questo motivo, il filtro viene spesso indicato come filtro di Wiener–Kolmogorov.

Il filtro di Wiener rappresenta il primo esempio di progettazione di un dispositivo di telecomunicazioni su base statistica; a questa scelta, e a quella quasi contemporanea di Shannon, che dette un analogo fondamento alla sua teoria dell'informazione, si deve la nascita della moderna teoria delle comunicazione come scienza statistica.

Dopo il filtro di Wiener sono stati elaborate molte tecniche statistiche di elaborazione dei segnali, tra i quali il filtro di Kalman.

L'ingresso di un filtro di Wiener si assume essere un segnale ${\displaystyle s(t)}$ alterato da un rumore additivo ${\displaystyle n(t)}$. L'uscita ${\displaystyle {\hat {s}}(t)}$ è calcolata per mezzo del filtro ${\displaystyle g(t)}$ usando la seguente convoluzione:

${\displaystyle {\hat {s}}(t)=g(t)*(s(t)+n(t))}$

dove il simbolo " * " indica l'operazione di convoluzione e:

• ${\displaystyle s(t)}$ è il segnale originale da ricostruire il più fedelmente possibile in uscita
• ${\displaystyle n(t)}$ è il rumore
• ${\displaystyle {\hat {s}}(t)}$ è il segnale stimato che speriamo sia uguale a ${\displaystyle s(t)}$
• ${\displaystyle g(t)}$ è il filtro di Wiener

L'errore è ${\displaystyle e(t)=s(t+\alpha )-{\hat {s}}(t)}$ mentre l'errore quadratico è ${\displaystyle e^{2}(t)=s^{2}(t+\alpha )-2s(t+\alpha ){\hat {s}}(t)+{\hat {s}}^{2}(t)}$ dove

• ${\displaystyle s(t+\alpha )}$ è l'uscita desiderata del filtro
• ${\displaystyle e(t)}$ è l'errore

A seconda del valore di α si affronta un tipo di problema diverso:

• Se ${\displaystyle \alpha >0}$ allora il problema è detto di predizione
• Se ${\displaystyle \alpha =0}$ allora il problema è detto di filtraggio
• Se ${\displaystyle \alpha <0}$ allora il problema è detto di smoothing

Introducendo esplicitamente l'espressione dell'integrale di convoluzione, ${\displaystyle {\hat {s}}(t)}$ può essere scritto come: ${\displaystyle {\hat {s}}(t)=g(t)*(s(t)+n(t))=\int _{-\infty }^{\infty }{g(\tau )\left[s(t-\tau )+n(t-\tau )\right]d\tau }}$.

e pertanto il valore atteso dell'errore quadratico sarà

${\displaystyle E(e^{2})=R_{s}(0)-2\int _{-\infty }^{\infty }{g(\tau )R_{x\,s}(\tau +\alpha )d\tau }+\int _{-\infty }^{\infty }{\int _{-\infty }^{\infty }{g(\tau )g(\theta )R_{x}(\tau -\theta )d\tau }d\theta }}$

dove

• ${\displaystyle \,\!x(t)=s(t)+n(t)}$
• ${\displaystyle \,\!R_{s}}$ è la funzione di autocorrelazione di ${\displaystyle s(t)}$
• ${\displaystyle \,\!R_{x}}$ è la funzione di autocorrelazione di ${\displaystyle x(t)}$
• ${\displaystyle R_{x\,s}}$ è la funzione di correlazione mutua, o cross-correlazione di ${\displaystyle x(t)}$ e ${\displaystyle s(t)}$

La ${\displaystyle g(t)}$ ottima è quindi quella che minimizza ${\displaystyle E(e^{2})}$.

Se il segnale ${\displaystyle s(t)}$ e il rumore ${\displaystyle n(t)}$ sono non correlati, ovvero la cross-correlazione è zero, l'espressione di semplifica in quanto

• ${\displaystyle R_{x\,s}=R_{s}}$
• ${\displaystyle \,\!R_{x}=R_{s}+R_{n}}$

Il filtro di Wiener ammette soluzione per due possibili casi: il caso in cui si voglia un filtro causale, e il caso in cui un filtro non causale sia accettabile. L'ultimo è più semplice, ma non si presta per applicazioni in tempo reale. L'obiettivo principale di Wiener era risolvere il caso in cui la condizione di causalità fosse valida.

${\displaystyle G(s)={\frac {S_{x,s}(s)e^{\alpha s}}{S_{x}(s)}}}$

Posto che ${\displaystyle g(t)}$ sia ottimale, allora l'equazione di minimo errore quadratico medio si riduce in ${\displaystyle E(e^{2})=R_{s}(0)-\int _{-\infty }^{\infty }{g(\tau )R_{x,s}(\tau +\alpha )d\tau }}$

e la soluzione ${\displaystyle g(t)}$ è l'inversa della trasformata di Laplace a due lati di ${\displaystyle G(s)}$.

${\displaystyle G(s)={\frac {H(s)}{S_{x}^{+}(s)}}}$

Dove

• ${\displaystyle H(s)}$ è la parte causale di ${\displaystyle {\frac {S_{x,s}(s)e^{\alpha s}}{S_{x}^{-}(s)}}}$ (cioè la parte di questa frazione avente una soluzione temporale positiva per trasformazione inversa di Laplace)
• ${\displaystyle S_{x}^{+}(s)}$ è la componente causale di ${\displaystyle S_{x}(s)}$ (p.e. la trasformata inversa di Laplace di ${\displaystyle S_{x}^{+}(s)}$ è non nulla solo se ${\displaystyle t\geq 0}$)
• ${\displaystyle S_{x}^{-}(s)}$ è la componente anti causale di ${\displaystyle S_{x}(s)}$ (p.e. la trasformata inversa di Laplace ${\displaystyle S_{x}^{-}(s)}$ è non nulla solo per t non negativo)

Questa espressione generale è complicata e necessita di una più dettagliata spiegazione. Per ottenere la soluzione ${\displaystyle G(s)}$ per un caso specifico, bisogna seguire i seguenti passi:^[6]

1. Si parte dallo spettro ${\displaystyle S_{x}(s)}$ in forma razionale e lo si fattorizza nelle componenti causale e anti causale:

${\displaystyle S_{x}(s)=S_{x}^{+}(s)S_{x}^{-}(s)}$

dove ${\displaystyle S^{+}}$ contiene tutti gli zeri e i poli nel semipiano sinistro (LHP) e ${\displaystyle S^{-}}$ contiene gli zeri e i poli nel semipiano destro (RHP). Si intendono semipiano destro e sinistro i due semipiani a sinistra e a destra dell'origine, essendo la retta verticale s=0 la frontiera dei due semipiani.

2. Si divide ${\displaystyle S_{x,s}(s)e^{\alpha s}}$ per ${\displaystyle {S_{x}^{-}(s)}}$ si espande il risultato in frazioni parziali.

3. Si selezionano quei termini dell'espansione che hanno poli nel semipiano LHP, termini indicati con ${\displaystyle H(s)}$.

4. Si divide ${\displaystyle H(s)}$ per ${\displaystyle S_{x}^{+}(s)}$. Il risultato è la funzione filtro di trasferimento desiderata ${\displaystyle G(s)}$

Per ricavare i coefficienti del filtro di Wiener si considera un signale w[n] passato a un filtro di Wiener di ordine N con coefficienti ${\displaystyle a_{i},i=0,\ldots ,N}$. L'uscita del filtro, indicata con x[n], è data dall'espressione

${\displaystyle x[n]=\sum _{i=0}^{N}a_{i}w[n-i]}$

L'errore residuo è indicato con e[n] e definito come e[n]=x[n]-s[n] (Vedi il corrispondente diagramma a blocchi). Il filtro di Wiener è concepito in maniera tale da minimizzare la deviazione standard (criterio MMSE), scritta concisamente nella maniera seguente:

${\displaystyle a_{i}=arg.min~E\{e^{2}[n]\}}$

dove ${\displaystyle E\{.\}}$ indica l'operatore di valore atteso. In generale i coefficienti ${\displaystyle a_{i}}$ sono essere complessi e possono essere derivati anche nel caso in cui w[n] e s[n] sono complesse. Per semplicità si considera solo il caso in cui tutte queste quantità sono reali. La deviazione standard può essere riscritta come:

${\displaystyle {\begin{array}{rcl}E\{e^{2}[n]\}&=&E\{(x[n]-s[n])^{2}\}\\&=&E\{x^{2}[n]\}+E\{s^{2}[n]\}-2E\{x[n]s[n]\}\\&=&E\{{\big (}\sum _{i=0}^{N}a_{i}w[n-i]{\big )}^{2}\}+E\{w^{2}[n]\}-2E\{\sum _{i=0}^{N}a_{i}w[n-i]s[n]\}\end{array}}}$

Per trovare il vettore ${\displaystyle [a_{0},\ldots ,a_{N}]}$ che minimizza l'espressione sopra si calcola la sua derivata rispetto ai coefficienti ${\displaystyle a_{i}}$

${\displaystyle {\begin{array}{rcl}{\frac {\partial }{\partial a_{i}}}E\{e^{2}[n]\}&=&2E\{{\big (}\sum _{j=0}^{N}a_{j}w[n-j]{\big )}w[n-i]\}-2E\{s[n]w[n-i]\}\quad i=0,\ldots ,N\\&=&2\sum _{j=0}^{N}E\{w[n-j]w[n-i]\}a_{j}-2E\{w[n-i]s[n]\}\end{array}}}$

Se si suppongono w[n] e s[n] stazionarie, si possono introdurre le seguenti sequenze ${\displaystyle R_{w}[m]~{\textit {and}}~R_{ws}[m]}$, note rispettivamente come l'autocorrelazione di w[n] e la cross-correlazione tra w[n] e s[n], definite come segue

${\displaystyle R_{w}[m]=E\{w[n]w[n+m]\}\quad {\textit {and}}\quad R_{ws}[m]=E\{w[n]s[n+m]\}}$

La derivata della MSE può dunque essere riscritta come (si noti che ${\displaystyle R_{ws}[-i]=R_{sw}[i]}$)

${\displaystyle {\frac {\partial }{\partial a_{i}}}E\{e^{2}[n]\}=2\sum _{j=0}^{N}R_{w}[j-i]a_{j}-2R_{sw}[i]\quad i=0,\ldots ,N}$

Imponendo la derivata uguale a zero si ottiene

${\displaystyle \sum _{j=0}^{N}R_{w}[j-i]a_{j}=R_{sw}[i]\quad i=0,\ldots ,N}$

che può essere riscritta in forma matriciale

${\displaystyle {\begin{bmatrix}R_{w}[0]&R_{w}[1]&\cdots &R_{w}[N]\\R_{w}[1]&R_{w}[0]&\cdots &R_{w}[N-1]\\\vdots &\vdots &\ddots &\vdots \\R_{w}[N]&R_{w}[N-1]&\cdots &R_{w}[0]\end{bmatrix}}{\begin{bmatrix}a_{0}\\a_{1}\\\vdots \\a_{N}\end{bmatrix}}={\begin{bmatrix}R_{sw}[0]\\R_{sw}[1]\\\vdots \\R_{sw}[N]\end{bmatrix}}}$

Queste equazioni sono note come equazioni di Wiener-Hopf. La matrice che appare nell'equazione è una matrice di Toeplitz simmetrica. Queste matrici sono definite positive e dunque non singolari, per cui esiste una sola soluzione nella determinazione dei coefficienti del filtro di Wiener. Inoltre esiste un algoritmo efficiente per risolvere le equazioni di Wiener-Hopf noto come algoritmo di Levinson-Durbin^[7].

Il filtro di Wiener è molto utile anche nel campo del riconoscimento di segnale, in particolare del Pile Up (sovrapposizione di due segnali nella stessa finestra temporale). Ecco come si esprime la funzione di trasferimento del filtro.

Si supponga di avere un sistema con una risposta r(t) al segnale unitario e un segnale non corrotto u(t) in ingresso (es: ${\displaystyle \delta (t-t_{0})}$). L'uscita s(t) del sistema sarà pertanto una convoluzione tra la funzione di risposta e il segnale di ingresso u(t):

${\displaystyle s(t)=\int _{-\infty }^{+\infty }{r(\tau )u(t-\tau )d\tau }}$

che nel dominio di Fourier si riduce semplicemente a ${\displaystyle s(\omega )=r(\omega )u(\omega )}$.

Nel caso reale invece l'uscita del sistema sarà un segnale corrotto c(t) = s(t) + n(t). La ${\displaystyle u(\omega )}$ in questo caso sarà:

${\displaystyle u(\omega )={\frac {c(\omega )}{r(\omega )}}}$

Lo scopo del Filtro di Wiener è quello di trovare un filtro ${\displaystyle H(\omega )}$ tale per cui vale:

${\displaystyle u'(\omega )={\frac {c(\omega )H(\omega )}{r(\omega )}}}$

dove ${\displaystyle u'(\omega )}$ è la miglior stima di ${\displaystyle u(\omega )}$, ovvero ${\displaystyle H(\omega )}$ minimizza la differenza quadratica tra U e U'. In formule equivale a minimizzare

${\displaystyle L(\omega )=\int _{-\infty }^{+\infty }{\left|u'(\omega )-u(\omega )\right|^{2}d\omega }=\int _{-\infty }^{+\infty }{\left|{\frac {c(\omega )H(\omega )}{r(\omega )}}-{\frac {s(\omega )}{r(\omega )}}\right|^{2}d\omega }}$

Risolvendo ${\displaystyle dL/dH=0}$ si ottiene l'espressione cercata:

${\displaystyle H(\omega )={\frac {|s(\omega )|^{2}}{|s(\omega )|^{2}+|n(\omega )|^{2}}}}$

Se si considera la u(t) come una delta, l'espressione si semplifica:

${\displaystyle H(\omega )={\frac {r^{*}(\omega )}{|r(\omega )|^{2}+|n(\omega )|^{2}}}\qquad u'(\omega )={\frac {c(\omega )r^{*}(\omega )}{|r(\omega )|^{2}+|n(\omega )|^{2}}}}$

Questo filtro può tornare molto utile nella ricerca di segnali, come per esempio in un algoritmo di trigger o riconoscimento di pattern, dato che la U' estrapolata rappresenta appunto una delta, quindi è facilmente riconoscibile.

• (EN) Robert Grover Brown e Patrick Y. C. Hwang, Introduction to Random Signals and Applied Kalman Filtering, 3ª ed., New York, John Wiley & Sons, 1996, ISBN 0-471-12839-2.
• (EN) Thomas Kailath, A View of Three Decades of Linear Filtering Theory, in IEEE Trans. Inf. Theory, vol. 20, 1974, pp. 146-181, DOI:10.1109/TIT.1974.1055174.
• (RU) A. N. Kolmogorov, Sequenze stazionarie in spazi di Hilbert, in Boll. Univ. Mosca, vol. 2, n. 6, 1941, pp. 1-40.
• (EN) N. Wiener, Extrapolation, Interpolation, and Smoothing of Stationary Time Series, New York, John Wiley & Sons, 1949, ISBN 0-262-73005-7.

• Norbert Wiener
• Filtro di Kalman

• Wikimedia Commons contiene immagini o altri file su Filtro di Wiener

Portale Matematica: accedi alle voci di Wikipedia che trattano di matematica