Formalismo ed applicabilità
del metodo ICA
(Independent Component Analysis)
Francesca Marcucci
Università di Perugia e INFN
Udine 31 gennaio 2003
ICA: Independent Component Analysis
ICA è una tecnica statistica per la decomposizione di un complesso dataset
nelle sue sottoparti indipendenti ed è particolarmente utile nella soluzione
di problemi di Blind Source Separation (BSS)
Modello:
Supponiamo che  M segnali di media nulla s1 s2 s3 ....... sM ma siano osservabili solo
N combinazioni lineari delle variabili si
xj = aij
· si
i=1,..,M
j=1,..,N
x=As
A={aij}
Statistical “latent variables” model
Se la matrice A non è nota il problema puo’ essere risolto facendo alcune assunzioni
sulle proprietà statistiche delle sorgenti si
Dovrebbe essere considerato anche un termine aggiuntivo n per il rumore
x=As + n
Ipotesi per l’ applicabilità di ICA:
 N >= M
(nel seguito assumiamo M=N senza perdita di generalità)
 al massimo una delle sorgenti e’ gaussiana
 le sorgenti si sono statisticamente indipendenti
 la matrice A ha rango massimo
 per ora n=0 , ma il modello puo’ essere esteso anche se di più difficile
risoluzione
Ambiguità del metodo:
Il metodo fornisce una misura dell’ indipendenza delle componenti ma non da
informazioni sull’Energia (varianza) e sull’ ordine in cui si ottengono, ovvero
la matrice A puo’ essere scritta (dopo la convergenza) come:
A=PD
Soluzione: Whitening or Sphering
P=permutazione D=matrice diagonale
PREPROCESSAMENTO DEI DATI:
Se i dati si non hanno media nulla allora si sottrae il valor medio
(ad xi)
Whitening o sphering:
Serve ad ottenere dei nuovi dati x con varianza unitaria
x=Vx dove E{ x xT }=1 (=I)
Se E{ x xT }=C allora V=C-1/2
infatti E{ x xT }= E{V x xT VT}=C-1/2 C C-1/2=I
Illustrazione del metodo:
Supponiamo di avere due variabili indipendenti uniformemente distribuite nella
regione illustrata, con media nulla e varianza unitaria
1/23
|si|<3
s2
Ad es. P(si)=
0 altrove
s1
Applichiamo
A=
2
3
2
1
Le direzioni
x2
ci danno
informazione sulle colonne di A
x1
… cerchiamo un metodo piu’ generale
Procedimento:
Si basa sul teorema del “limite centrale”:
La distribuzione della somma di variabili random indipendenti tende ad
una distribuzione gaussiana
Per stimare una delle componenti indipendenti consideriamo
y = wT x = i wi xi
se w fosse l’ i-ima riga di A-1 allora y= si
z = AT w

y = wT x = wT As = zT s
Se i dati hanno varianza unitaria W e’ una matrice ortogonale WWT=I
Come usare il teorema del limite centrale?
Ora abbiamo y = zTs ossia una combinazione lineare delle
sorgenti indipendenti. Tale somma è “piu’ gaussiana” delle
componenti originarie e lo diventa “al minimo” quando y=si ossia z ha solo
l’i-imo elemento non nullo.
w scelto in modo da massimizzare la non-gaussianità di wTx
Misure di non-gaussianità:
KURTOSIS
kurt(y)=E{y4} – 3(E{y2})2
è nullo per variabili gaussiane
quindi si cerca il max di |Kurt(y)|
NEGENTROPY J(y)=H(ygauss) – H(y)
con H(y)= f(y) log f(y) dy
è nulla per variabili gaussiane (quelle con la max entropia H)
MUTUAL INFORMATION
I(y1,…,yM) = i H(yi) – H(y)
È nulla per variabili indipendenti e non negativa  va minimizzata
Modello di rete neurale:
W
Q

x
x
y
y è una stima del vettore s
y=Wx
Q è una stima della matrice A
x=Qy

1. “Apprende” una matrice W tale che y=Wx sono indipendenti
2. “Apprende” una matrice Q tale da minimizzare
E{||n’||2}=E{||x-Qy||2}
W=BTV
Con pre-withening:
V
Q
BT
x

x
x
y
LEARNING:
Massimizzare/minimizzare rispetto a w una delle funzioni F(w) precedenti
imponendo dei vincoli ad esempio E{y2}=1 e E{y}=0 , ad esempio utilizzando i
moltiplicatori di Lagrange :
gradient-ascendent method:
wk+1 = wk +  L ’wk
Newton-Like method:
L ’’w = r(w) Rxx
2
L ’’w wk= -L ’wk
2
wk+1 = wk -  Rxx-1 L ’wk/ r(w)
ALGORITMI:
Herault-Jutten: fallisce per piu’ di 2 sorgenti
EASI: performance uniforme
Bell’s and Seinowsky’s: performance uniforme e non richiede pre-withening
Chicocki and Amari: per feedforward e recurrent network
BIGRADIENT: necessario prewhitening, molto flessibile
NONLINEAR PCA: senza prewhitening separa solo componenti sinusoidali.
Adatto principalmente per funzioni sub-gaussiane
FastICA :
Caso semplice “one-unit”
(una sola unità computazionale  1 neurone con peso w)
FastICA trova un vettore unitario w tale che massimizzi la non-gaussianità
di wTx (utilizzando la Negentropy) con il metodo Newton-Like
1. Sceglie un iniziale vettore w random
2. calcola w+ = E{xg (wTx)} – E{g ’(wTx)} w
g derivata di una funzione
non quadratica
3. controlla se w = w+ / ||w+||
4. se non converge (w  w+ = 1, hanno la stessa direzione) ritorna al punto 2
Tale algoritmo one-unit permette di determinare solo 1 componente ma può
essere facilmente esteso per la stima di più componenti indipendenti
improntando una rete “several-unit” con neuroni di pesi w1,…,wn
Converge più rapidamente del metodo ICA; non necessita della stima di funzioni g o
di parametri di altri parametri, è gratuito e disponibile sul web.
Recente applicazione:
(Baccigalupi et al. 2002)
SCOPO: Separazione di componenti astrofisiche sovrapposte,
ricostruendone sia le caratteristiche spaziali che spettrali, senza
assunzioni a priori se non l’indipendenza e l’assenza di componenti gaussiane
MODELLO:
xi(r,)=ij sj(r,)
(N differenti processi fisici)
x=vettore M-dim
M canali di misura (diverse bande di frequenza) e strumento caratterizzato
da una PSF B(r,) e funzione di risposta t(’)
x(r)= B(r-r’,’)j t(’)sj(r’,’) dr’ d’ + n
Ipotesi:
1. funzione separabile sj(r,) = fj() sj(r)
2. B(r-r’,)=B(r) indipendente dalla frequenza
3. aij= t(’)f (’) d’

x(r)=A s(r) * B(r) + n
4. n è un rumore bianco, indipendente dal segnale , Gaussiano e stazionario
Synchrotron
angular
power
spectra
Input
output
CMB
angular
power
spectra
Input
output
Limiti:
La ricostruzione della matrice di separazione
peggiora nellì ipotesi in cui il rapporto tra due
componenti è fortemente variabile lungo la skymap
ES: le polveri dominano sul piano galattico mentre
CMB domina ad alte latitudini
La ricostruzione è ottenuta con un’ errore migliore
dell’ 1% nelle regioni in cui S/N  1.5 ,
l’errore cresce fino al 10% per S/N  1
Ancora un’ applicazione in astrofisica:
(Maria Funaro, Erkki Oja e Harri Valpola,2002)
Scopo: Individuare e rimuovere gli “artefacts” che influenzano le immagini
(fluttuazioni,stelle della nostra galassia, rumore strumentale) basandosi
sull’ analisi del profilo temporale della luminosità dei pixel e
sull’indipendenza delle componenti dell’ immagine.
Dati: Immagini della Galassia M31
Modello:
N pixel
X matrice TxN
T immagini
M sorgenti
X = AS
riga Xt: singola immagine al tempo t
colonna Xn: serie temporali (curve di luce) del pixel n
 S matrice MxN
righe: immagini delle componenti indipendenti
per il singolo pixel n Xn =
 am Smn
 A matrice TxM : le M colonne di A (mixing vectors am) sono delle “curve di luce
virtuali” le cui combinazioni lineari danno quelle reali Xn
am caratterizza il comportamento temporale della sorgente m
Smn caratterizza il comportamento spaziale
T = 35
e N=100x100 pixel
dopo whitening 3510 componenti indip.
Immagine originaria
1° e 2° autovettori: Raggi cosmici
5° autovettore: Sorgente puntiforme
Conclusioni

Ci sono pochi casi in letteratura di applicazioni ICA in
astrofisica, ma in questo campo l’indipendenza delle componenti
assicura l’applicabilità del metodo.
 La bontà statistica del metodo e’ legata principalmente alla
minimizzazione della funzione di costo nella rete neurale
implementata .
E’ necessario verificarne l’accuratezza con modelli simulati piu’
vicini alla realtà osservativa
 ICA è sicuramente piu’ rapido dei metodi tradizionali … ma è
ugualmente attendibile?
PROPOSTA:
 Utilizzare in un primo momento FastICA ,Likelihood e Wavelet
con gli output del light simulator e confrontare i due metodi