Tema 2: Istogrammi, trasformazioni di v.a.,
e misura di indici statistici
Esempio di v.a.: Immagini radar
Baltrum, isole Frisone
Immagine grezza
(con disturbo aleatorio)
Immagine filtrata
Concetto di Istogramma
x
x 

Pi  Pr  xi 
 X  xi  
2
2 

xi x 2

xi x 2
K ( xi )
f X ( x)dx  f X ( xi )x 
N
K(xi) = Num. osservati che cadono
nell’intervallo (bin) centrato in xi
Istogramma ottenuto mediante
N=1000 realizzazioni,
{ x[k]; k=1,2,..,N }, di una v.a.
Gaussiana standard
200
180
160
Normalizzazione:
140
120
1 K ( xi )
f X ( xi )  
N x
100
80
60
40
Nota sulla accuratezza della misura:
deve essere soddisfatta la regola empirica
N  100 / Pi dove Pi  f X ( xi )x
20
0
-4
-3
-2
-1
0
1
2
3
4
“bin”
x  largh. bin
Istogramma di immagini radar
Istogramma non normalizzato
Dettaglio immagine grezza
900
800
700
delle ampiezze
La mappa riporta il logaritmo
1000
v.a. ampiezza
600
500
400
300
200
100
0
0
1000
3000
4000
5000
6000
Istogramma normalizzato
Istogramma normalizzato
-7
8
-4
7
2000
x 10
x 10
Trasformaz.
v.a. ampiezza
6
7
v.a. intensità
6
5
pA (a) 
4
a
2
e

a2
2 2
di v.a.:
u (a )
IA
2
5
p I (i ) 
4
1
2
e

i
2
u (i )
3
3
2
2
1
1
0
0
0.5
1
1.5
2
2.5
3
3.5
7
0
x 10
0
1000
2000
3000
4000
5000
6000
Nota: disturbi aleatori con ddp analoghe si trovano anche nei ponti radio
Calcolo di istogrammi
- Esercizio: Calcolare e visualizzare l’istogramma normalizzato
utilizzando N=105 realizzazioni di una v.a. Gaussiana standard,
e confrontarlo con la ddp ideale:
 x2 
1
f X ( x) 
 exp   
2
 2
- Problema: fissato N, come scegliere x (ovvero il numero di bin)?
Verificare il compromesso risoluzione/accuratezza
1 K ( xi )
f X ( xi )  
N x
[ Istruzioni utili: randn, function, size, hist, bar ]
Esempio di file.m: istog.m
% Calcolo istogramma normalizzato
function [ddp,x0] = istog(NBIN,x)
% IN:
% OUT:
%
%
vettore dati, x; numero di bin, nbin;
istogramma normalizzato, ddp;
vettore valori di centro-bin, x0;
uscita su video di istogramma normalizzato
Nx=size(x,2);
% numero di dati
[count,x0]= hist(x,NBIN);
Deltax=x0(2)-x0(1);
ddp= count/Deltax/Nx;
% esegue i conteggi
% calcola larghezza bin
% normalizza
bar(x0,ddp)
% grafico
Esempio di risultati
% Compromesso accuratezza/risoluzione dell'istogramma
x=randn(1,10000);
figure
[d,x0]=istog(10,x);
figure
[d,x0]=istog(50,x);
figure
[d,x0]=istog(500,x);
0.5
0.45
0.4
0.35
0.3
0.25
0.2
0.45
0.15
0.4
0.1
0.35
0.05
0.3
0
-4
0.4
-3
-2
-1
0
1
2
3
0.25
0.35
0.2
0.3
500 bin
0.15
0.25
0.1
0.2
0.05
0
-4
0.15
-3
-2
-1
0
1
2
3
4
0.1
50 bin
0.05
0
-4
-3
-2
-1
0
10 bin
1
2
3
4
[ per confronto si è sovrapposta la ddp effettiva ]
4
Generazione di v.a. exp. negativa
Generazione di v.a. Y con d.d.p. fY(y) e funzione di distribuzione FY (y):
Y  FY1( X ) con X v.a. uniforme su [0,1]
- Generare N=104 realizzazioni di una v.a. esponenziale negativa
con parametro h=1
[ Istruzioni utili: rand, ... ]
- Calcolare e visualizzare l’istogramma normalizzato dei dati generati
[ Istruzioni utili: hist, bar ]
Esempio di file.m: genespneg.m
% Generazione v.a. exp. negative con parametro eta=1
function X = genespneg(N)
% IN: numero di v.a., N;
% OUT: vettore di v.a. esponenziali negative (eta=1);
unif=rand(1,N);
% genera N v.a. uniformi tra 0 e 1
X=-log(1-unif);
%
%
%
%
trasformazione di v.a. per ottenere
le v.a. esponenziali negative (eta=1);
nota: si può usare anche la trasformazione
1
X=-log(unif)
0.9
0.8
0.7
0.6
% generazione 20000 v.a.
x= genespneg(20000);
% istogramma normalizz. con 100 bin
[d,x0]=istog(100,x);
0.5
0.4
0.3
0.2
0.1
0
0
2
4
6
8
10
12
Misura del valor medio e deviazione standard
Dati: { x[k]; k=1,2,..,N }
Nbin

h X  E{ X }   xf X ( x)dx   xi f X ( xi )x

i 1
Nbin
xi K ( xi ) 1
N

  k 1 x[k ]  hˆX
N
N
i 1
 X  E{( X h X ) }  ˆ X 
2
K ( xi )
Pi  f X ( xi )x 
N
1
N
2
ˆ
 x[k ] h X 

k 1
N
- Valutare sperimentalmente il valor medio e la deviazione standard
della v.a. esponenziale negativa generata,
e di una v.a. uniforme nell’intervallo [-1,3]
di cui si hanno N=1000 e N=10 realizzazioni generate al calcolatore
[ Istruzioni utili: rand, mean, std ]
Esempio di file.m & risultati (1)
% Generazione delle v.a. esponenziali negative con eta=1
x=genespneg(20000);
% calcola media e deviazione standard
mean(x)
std(x) % oppure std(x,1)
ans =
1.0102
ans =
1.0034
Valori effettivi: media=1; dev_stand=1
Esempio di file.m & risultati (2)
% Generazione 1000 realizzazioni di v.a. uniforme tra –1 e 3
x=rand(1,1000)*4-1;
% calcola media e deviazione standard
mean(x)
std(x) % oppure std(x,1)
% calcolo con sole 10 realizzazioni
mean(x(1:10))
std(x(1:10)) % oppure std(x(1:10),1)
ans =
1.0645
ans =
1.1713
N=1000
ans =
1.2035
ans =
0.8916
Valori effettivi: media=1; dev_stand=1.1547
N=10
Skewness e Kurtosis
Definizione di Skewness e Kurtosis:
3 
 X (3)
32
  X (2) 
4 
 X (4)
3
2
  X (2) 
dove  X (2)   X2 è la varianza
Sono indici di non Gaussianità !
1
N
n
ˆ
 X (n)  E{( X  h X ) }  ˆ X (n)   k 1  x[k ]  h X 
N
n
Misura di Indici Statistici
Esercizio:
• Generare una sequenza di N=105 realizzazioni di una v.a. X, Gaussiana
standard e di una v.a. Y, exp monolatera, avente varianza unitaria
• Misurare dai dati valor medio, varianza e skewness al variare del numero n
di campioni utilizzati per la misura, e graficare il risultato per 1<n<N
Osservazione: si noti il diverso comportamento al variare di n, in particolare si
noti come all’aumentare dell’ordine del momento misurato il numero di dati
necessario per ottenere una misura “accurata” aumenti …
Esempio di risultati (3)
Misure empiriche (stime) di valor medio, dev. standard, skewness e kurtosis
al variare del numero n di campioni utilizzati per la misura: n=1,2, …, N
3  0
4  0
Esempio di risultati (4)
3  2
4  6
Esempio di file.m: skewkurt.m
% Misura della media, dev. standard, skewness e kurtosis di v.a. Gaussiane ed esponenziali
clc
clear
N=10000;
EN=exprnd(1,1,N); % genera N v.a. esponenziali negative (eta=1);
GA=randn(1,N);
% genera N v.a. Gaussiane standard;
for indice=1:N/10,
ncamp=indice*10;
etaEN(indice)=mean(EN(1:ncamp));
etaGA(indice)=mean(GA(1:ncamp));
stdEN(indice)=std(EN(1:ncamp),1);
stdGA(indice)=std(GA(1:ncamp),1);
skeEN(indice)=skewness(EN(1:ncamp));
skeGA(indice)=skewness(GA(1:ncamp));
kurEN(indice)=kurtosis(EN(1:ncamp))-3;
kurGA(indice)=kurtosis(GA(1:ncamp))-3;
% per numero di campioni variabile a passi di 10
% valuta le medie
% valuta le deviazioni standard
% valuta le skewness
% valuta le kurtosis
end
a=10*[1:N/10];
% esegue i grafici
figure
semilogx(a,etaEN,'r',a,stdEN,'g',a,skeEN,'b',a,kurEN,'k')
xlabel('num. campioni')
title('v.a. esponenziale negativa (eta=1)')
legend('media','dev. stand.','skew.','kurt.')
figure
semilogx(a,etaGA,'r',a,stdGA,'g',a,skeGA,'b',a,kurGA,'k')
xlabel('num. campioni')
title('v.a. Gaussiana standard')
legend('media','dev. stand.','skew.','kurt.')