Features extraction nel dominio della frequenza Fourier Transform (FT)

1. Features extraction nel dominio della frequenza • Fourier Transform (FT) • Short-Time Fourier Transform (STFT) • Continuous Wavelet Transform (CWT) • Discrete Wavelet Transform (DWT) • Matching Pursuit (MP)

2. Fourier Transform Durante il 19-esimo secolo Fourier dimostrò che ogni funzione periodica può essere espressa come somma di funzioni esponenziali complesse: Il segnale x(t) è moltiplicato per un termine esponenziale complesso caratterizzato da una particolare frequenza f, quindi integrato su tutti gli istanti di tempo. Se il segnale presenta componenti di frequenza f, tale componente “matcherà” con il termine sinusoidale ed il prodotto tra di essi avrà valore (relativamente) alto.

3. Segnale stazionario FT Segnale non stazionario FT

4. Fourier Transform FT non è in grado di distinguere i due segnali in maniera efficiente (lo spettro è molto simile, nonostante i segnali siano profondamente differenti!). Per FT i segnali sono molto simili dal momento che essi presentano le medesime componenti frequenziali. FT evidenzia la presenza delle componenti armoniche, ma non permette di avere informazioni su dove tali frequenze siano posizionate (temporalmente). Pertanto, FT non è uno strumento valido nello studio di segnali non-stazionari (la maggior parte dei segnali reali fisici), ovvero segnali le cui componenti cambiano nel tempo. Se non siamo interessati agli istanti di tempo in cui le componenti frequenziali si presentano, ovvero siamo interessati a conoscere se tali componenti sono presenti (o meno) nel segnale, allora FT può tornare utile.

5. Considerazioni: La funzione di base e-j2πft è valutata lungo tutto l’asse temporale, per cui ci si limita ad un’analisi globale del segnale. D’altra parte, un’analisi locale del segnale permette una localizzazione nel tempo delle componenti spettrali. Per i segnali non stazionari occorre inserire nella trasformazione una dipendenza dal tempo. Il modo più immediato consiste nel rendere locale la FT, non operando su tutto il supporto del segnale x(t) ma su porzioni di esso, ottenute moltiplicando x(t) per una finestra che trasla nel tempo.

6. Short-Time Fourier Transform • Moltiplicando per g(t-τ) si limita lo spettro del segnale all'interno della finestra temporale definita dalla funzione g(t-τ) stessa. Si ottiene un'informazione sul suo contenuto armonico in un intorno dell'istante di tempo τ. • g*(t) è il complesso coniugato della funzione finestra (impulso gaussiano, impulso rettangolare, etc). • STFT mappa la funzione 1D x(t) in un dominio 2D (tempo/frequenza). • STFT fornisce lo spettro del segnale alterato dalla presenza della finestra (la finestra ha dimensioni prefissate -> la dimensione della finestra è un parametro critico!). • Vantaggi: # STFT fornisce una descrizione tempo/frequenza del segnale. # STFT supera le difficoltà della FT grazie al segnale finestra. Per ogni τ ed f è calcolato un coefficiente della STFT.

7. g*(t-τ2) g*(t-τ3) g*(t-τ1) τ1 τ2 τ3 La funzione Gaussiana rossa individua la finestra centrata in t=τ1, la blu quella in t=τ2, la verde in t=τ3. A ciascuna finestra corrisponde una STFT (in 3 istanti di tempo differenti). In questo modo otteniamo una vera rappresentazione tempo/frequenza del segnale.

8. Nel segnale ci sono 4 componenti in frequenza. Il range 0-250ms presenta una sinusoide di 300 Hz, i successivi intervalli (tutti di 250ms) presentano sinusoidi di 200 Hz, 100 Hz e 50 Hz. Apparentemente è un segnale non-stazionario. Il grafico è simmetrico rispetto all’asse delle frequenze (la FT di un segnale reale è sempre simmetrica); dal momento che la STFT è la versione “finestrata” della FT, non ci sorprende che la STFT sia simmetrica rispetto all’asse delle frequenze. Notare la presenza dei 4 picchi corrispondenti alle 4 componenti in frequenza. Diversamente dalla FT, questi 4 picchi sono localizzati a differenti intervalli di tempo lungo l’asse del tempo (il segnale originale presenta 4 componenti spettrali localizzati in intervalli di tempo differenti). Contributo in bassa frequenza Contributo in alta frequenza

9. Grandi finestre: Buona risoluzione in frequenza (separazione delle armoniche), limitata risoluzione nel tempo. Finestre strette: Buona risoluzione nel tempo (capacità di localizzazione), limitata risoluzione in frequenza. Una finestra stretta è abile a rivelare componenti in alta frequenza, ma fornisce una risposta insoddisfacente per le componenti di bassa frequenza. Dal principio di indeterminazione di Heisenberg: se la finestra g(t) ha un supporto “largo”, allora la FT ha un supporto stretto (e viceversa) -> E’ impossibile ottenere contemporaneamente da una sola finestra una buona risoluzione nel dominio del tempo e nel dominio della frequenza! Quindi si possono conoscere gli intervalli di tempo in cui certe bande di frequenze esistono: Problema della Risoluzione STFT :: considerazioni

10. Finestra “piccola” g(t)=exp(-a(t2)/2) (a determina la lunghezza della finestra) Finestra “grande” Finestra “media”

11. Continuous Wavelet Transform Supera i problemi della STFT usando finestre di lunghezza variabile (finestre piccole per “catturare” il contenuto di alta frequenza garantendo buona localizzazione temporale; finestre grandi per “catturare” il contenuto in bassa frequenza garantendo buona risoluzione frequenziale) con x(t) = segnale = parametro di traslazione s = parametro di scaling = 1/f = wavelet madre se e solo se sono soddisfatti i seguenti vincoli:

12. Continuous Wavelet Transform Le funzioni kernel usate nella trasformata wavelet sono tutte ottenute da una funzione prototipo (wavelet madre), mediante operazione di scaling e/o traslazione Con a = parametro di scaling b = parametro di traslazione Continuous Wavelet Transform: Valori alti del parametro di scaling “a” permettono un’analisi grossolana (non di dettaglio) del segnale. Viceversa, bassi valori per “a” permettono di “catturare” l’informazione di dettaglio del segnale. D’altra parte, nel dominio di Fourier, a basse frequenze (alta scala) corrispondono informazioni globali, mentre ad alte frequenze (bassa scala) è associata l’informazione di dettaglio (o il rumore).

13. Analisi di dettaglio (bassa scala) In blu wavelet a supporto compatto con parametro di scaling unitario (bassa scala, alta frequenza). In figura la funzione wavelet è centrata in to=2, 40, 90 e 140. Per ogni versione traslata, la wavelet è moltiplicata per il segnale (in giallo): il prodotto è non nullo solo dove il segnale cade nella regione di supporto della wavelet; altrimenti è zero. Shiftando la wavelet il segnale è localizzato nel tempo; cambiando il parametro di scaling, il segnale è localizzato in frequenza (in scala).

14. Analisi grossolana (alta scala) Se il segnale ha componenti spettrali in corrispondenza del parametro di scaling, il prodotto della wavelet per il segnale è relativamente alto. Viceversa, se la componente spettrale corrispondente al parametro di scaling non è presente nel segnale, il valore del prodotto sarà relativamente basso.

15. Segnale composto da 4 componenti frequenziali a 30 Hz, 20 Hz, 10 Hz e 5 Hz. NB: Bassa scala -> alta frequenza (e viceversa) Informazione in bassa frequenza Informazione in alta frequenza

16. A piccoli valori del parametro di scaling corrispondono alte frequenze (la frequenza decresce all’aumentare del valore di scaling). La porzione di grafico con parametro di scaling vicino allo zero corrisponde alla più alta componente frequenziale, mentre ad alti valori di scala corrispondono contenuti di bassa frequenza. Osserviamo che il segnale presenta prima la componente di 30 Hz (la più alta frequenza) e ciò si manifesta a valori bassi di scaling e alla traslazione tra 0 e 30. Successivamente è presente il contenuto a 20 Hz, la seconda frequenza più alta, e cosi via. La componente a 5 Hz si manifesta verso la parte terminale dell’asse di traslazione, con alti valori di scala (bassa frequenza).

17. Inverse Wavelet Transform - Synthesis E’ possibile ricostruire il segnale originario partendo dai coefficienti wavelet con: Soggetta alla seguente condizione di ammissibilità della costante: FT di (t)

18. Esempi di Wavelets madri (t)

19. Analisi multirisoluzione - Discrete Wavelet Transform Il segnale originale x[n] subisce un filtraggio passa-alto (filtro g[n]) e passa-basso (filtro h[n]). L’output dei due filtraggi è sottocampionato (di 2). Tale attività completa la decomposizione wavelet di primo livello.

20. Analisi multirisoluzione - Discrete Wavelet Transform L’analisi Wavelet è un’estensione dell’analisi di Fourier. Si trasforma il segnale di input in un sottospazio rappresentato da versioni scalate e traslate della funzione di base denominata mother wavelet. S cD ~ 30 coeff. S 60 Samples cA1 cD1 cA ~ 30 coeff. cA2 cD2 Haar wavelet cA3 cD3

21. Matching Pursuit e dizionari overcompleti

22. Lo scopo dell’attività di rappresentazione dei dati (o segnali) è quello di individuare un appropriato insieme di coefficienti (features) capaci di descrivere numericamente (e possibilmente in maniera compatta) i dati di input Processo di classificazione CLUSTER Trasformazione COEFFICIENTI (FEATURES) SEGNALE Processo di sintesi SEGNALE SINTETIZZATO

23. Normalmente, nell’ambito delle attività legate al Pattern Recognition, sono usate PCA, Trasformata di Fourier e Trasformata Wavelet al fine di ottenere una rappresentazione del problema in questione Ora siamo interessati a studiare sistemi non-ortogonali di funzioni (over-completi, non basi) in cui la rappresentazione di ciascun elemento dello spazio di input non sia necessariamente unica. Si richiede che la rappresentazione possa essere ottenuta in modi diversi partendo da un insieme di funzioni linearmente dipendenti

24. Consideriamo un dizionario di poche migliaia di vocaboli per cui un particolare concetto può essere descritto usando i termini del vocabolario a fronte di ricorrere a frasi particolarmente lunghe. Viceversa, usando un dizionario di molti più vocaboli (ad esempio 100.000) lo stesso concetto può essere descritto brevemente, con frasi molto corte o, addirittura, con un solo vocabolo In generale, “dizionari” di grandi dimensioni (over-completi) permettono di descrivere lo stesso concetto in modi differenti per cui, tra le diverse rappresentazioni, è possibile individuare quella più opportuna per il particolare contesto applicativo. Si vuole individuare un insieme di funzioni che non siano necessariamente tra loro ortonormali o linearmente indipendenti, ma che possano costituire un operatore di rappresentazione usando espansioni di funzioni ridondanti

25. Rappresentazione di un segnale monodimensionale usando un insieme finito di funzioni linearmente dipendenti (ridondanti) {j} Il segnale f può essere recuperato usando l’equazione di sintesi Sia c il vettore ottenuto applicando l’operatore di frame F al segnale f (equazione di analisi) L’ultima relazione mette in evidenza la differenza tra l’operatore di frame F e l’operatore F†. F associa un vettore di coefficienti c (le features) al segnale decomponendolo (proiettandolo) in termini di atomi del dizionario. Viceversa, F† (operatore di sintesi o ricostruzione) permette il recupero del segnale come sovrapposizione degli atomi del frame duale opportunamente pesati a coefficienti c. Le colonne di F† sono gli elementi del frame duale.

26. Funzioni di Gabor monodimensionali La funzione di Gabor 1D è rappresentata da una funzione Gaussiana modulata da una funzione esponenziale complessa Esponenziale complesso Funzione Gaussiana La funzione di Gabor permette un’ottima localizzazione sia nel dominio del tempo che delle frequenze La funzione di Gabor può essere usata per rappresentare segnali e per estrarre features in un processo di classificazione

27. Dizionario overcompleto 1D di Gabor

28. Dizionario overcompleto 1D di Gabor

29. Dizionario overcompleto 1D di Gabor

30. Dizionario overcompleto 1D di Haar

31. 1 2 3 4 … … … … … L L’algoritmo di Matching Pursuit per l’estrazione di features MP decompone un segnale f in una espansione lineare di funzioni j che appartengono al dizionario overcompleto e seleziona gli atomi che meglio si adattano alla struttura globale di f MP è un algoritmo sub-ottimo che raffina l’approssimazione del segnale progressivamente mediante una procedura iterativa. MP ottimizza l’energia acquisita ad ogni iterazione Sia u0 = f ed L la cardinalità del dizionario F = {j} For k = 0 to n = numero di funzioni di Gabor desiderate Calcola per ciascun elemento di F: cj= <uk,j> , j = 1, ..., L Segnale uk (k = 0 alla prima iterazione) Seleziona l’atomo jk che massimizza jk =argmax cj j =1,…,L Inserisci jk nel sotto-dizionario F* Aggiorna il segnale residuo: uk+1 = uk - cjjk c1 = <uk, 1 > End loop c2 = <uk, 2 > c3 = <uk, 3 > c4 = <uk, 4 > Atomi del dizionario overcompleto … jk = argmax cj uk+1 = uk - cj jk … Selezione dell’atomo migliore Aggiornamento del segnale residuo … … … cL = <uk, L>

32. L’algoritmo di Matching Pursuit per l’estrazione di features Lo schema identifica solo un sotto-dizionario di pochi atomi in grado di rappresentare qualsiasi segnale monodimensionale Sotto-dizionario F* (Insieme di 256 funzioni del dizionario overcompleto) SCHEMA MP Rappresentazione di qualsiasi segnale mediante proiezione sul sotto-dizionario F* Si usa l’equazione di analisi della Teoria dei frame

33. L’algoritmo di Matching Pursuit con l’uso del dizionario di Gabor Errore di ricostruzione per il segnale usato nel processo di generazione del sotto-dizionario F* usando lo schema di MP (equazioni di analisi e sintesi) I primi atomi di F* catturano le informazioni di bassa frequenza del segnale I successivi atomi descrivono l’informazione di dettaglio (contenuto in alta frequenza) L’errore residuo converge a zero in un numero finito di iterazioni

34. L’algoritmo di Matching Pursuit con l’uso del dizionario di Gabor Errore di ricostruzione per un segnale non usato dallo schema MP per la generazione di F* (equazioni di analisi e sintesi) Si garantisce ancora una buona approssimazione del segnale Il sotto-dizionario F* spanna bene lo spazio dei segnali 1D Alta capacità di sintesi di F*

35. L’algoritmo di Matching Pursuit con l’uso del dizionario di Haar Errore di ricostruzione per il segnale usato nel processo di generazione del sotto-dizionario F* usando lo schema di MP (equazioni di analisi e sintesi) La convergenza dell’errore di ricostruzione è più veloce rispetto al caso di dizionari di Gabor

36. L’algoritmo di Matching Pursuit – Versione generalizzata Let R0p = fp with 1 p m and m=number of training signals Loop: for i = 0 to n=number of desired basis elements Compute Cg for each atom { i } in the dictionary Select the i-th basis element i that maximizes Cg Update coefficients for each block pi = < Rifp, i > Update the residue blocks Ri+1fp = Rifp - pi i Increment the iteration counter i = i + 1 End loop At iteration i let (p0, …, pi-1) the coefficients computed through iteration i-1 to represent block fp with 1 p m Compute the centroid  of (p0, …, pi-1, < Rifp, i >) Compute the mean distance U between  and (p0, …, pi-1, < Rifp, i >) Return U as value of Cg

39. RGB GREEN RED BLUE Esempio: An automated approach for vascular oral network detection Prevailing channel Few information Low contrast Best vessel/background contrast

40. GREEN RESTRICTED DYNAMICS HE POOR QUALITY CLAHE STRETCHED DYNAMICS

41. 2D Gaussian function Sinusoidal modulator 2D Gabor wavelet

42. 2D Gaussian function Sinusoidal modulator 2D Gabor wavelet

43. 2D Gaussian function Sinusoidal modulator 2D Gabor wavelet

44. CWT can be meant as correlation measure between the vessel pixels and a particular Gabor wavelet The intensity values of the equalized green channel are reversed so that blood vessels present values higher than the background pixels (otherwise CWT detects background!) CWT presents stronger responses when Gabor wavelet is tuned to vessels pixels (maximum correlation) A grey-levels image is defined by evaluating the scaling parameters for which the maximum modulus response is realized by the CWT In white low-frequency information (background) is shown. Darker regions are related to different size vessels (high frequency)

45. THRESHOLD Image binarization is realized according to the Otsu scheme[3] The automated threshold guarantees maximum distance of the clusters (background, vessel) and minimum intra-cluster variance [3] N. Otsu, "A Threshold Selection Method from Gray-Level Histograms," IEEE Transactions on Systems, Man and Cybernetics, vol.9, no.1, 1979, 62-66 VESSEL BK

46. Operatori morfologici Thinning Segmentazione finale

47. Esempio: Classificazione da immagini telerilevate mediante Gabor function