Tema 5: Misura di un segnale costante in rumore

1. Tema 5:Misura di un segnale costante in rumore Significato della ddp a priori versus la ddp a posteriori Modello scalare dell’osservato: campioni di rumore segnale costante aleatorio Modello vettoriale dell’osservato: Vettore valor medio degli osservati: Matrice di covarianza degli osservati: Varianza degli osservati:

2. Segnale costante in rumore: tre realizzzioni N.B. Vogliamo misurare il valore s assunto da S nella particolare realizzazione che viene osservata!! s=12.8 s=9.4 s=7.99

3. File.m: segncost_realiz.m % Grafico di una realizzazione del segnale osservato % “tempo-continuo” e tempo-discreto clcclear N=10; % numero di campioni prelevatisc=20; % fattore di sovracampionamento per simulare il "tempo continuo" etas=10 % valor medio del segnale costantesig2s=4 % varianza del segnale costantesig2w=4 % varianza del rumore s=randn(1,1)*sqrt(sig2s)+etas; % genera un valore di S con ddpGaussiana wtc=randn(N*sc,1)*sqrt(sig2w); % genera il rumore "tempo continuo" con ddp Gaussianaytc=s+wtc; % calcola il segnale "tempo continuo" totale y=ytc(1:sc:N*sc); % estrae gli N campioni osservati plot(ytc,'r') % traccia il grafico del segnale totale “tempo-continuo“xlabel('(msec)') % scala temporale – esemplificativa - %(tempo di campionamento= sc msec)ylabel('volt')axis([0 N*sc etas-5*sqrt(sig2s) etas+5*sqrt(sig2s)])hold on;plot([0 N*sc],[s s],'k') % traccia il grafico del livello costante;hold on;stem([1:sc:N*sc],y) % traccia i campioni osservati

4. Misura di un segnale costante in rumore Errore della misura Algoritmo di misura: Indici di qualità dell’algoritmo: valore medio e potenza media dell’errore Significato della ddp a priori: Ho a disposizione solo fS(s), ovvero nessun osservato; cosa possiamo fare? Significato della ddp a posteriori: Avendo a disposizione N osservati, cosa possiamo fare per migliorare l’accuratezza della misura?

5. Misura del segnale costante in rumore Algoritmo B: minimo errore quadratico medio Algoritmo A: media campionaria Potenza media statistica dell’errore di misura: Valore medio dell’errore di misura: Nota bene: entrambi sono combinazione lineare di v.a. Gaussiane, quindi entrambe sono v.a. Gaussiane

6. Algoritmo B: calcolo dell’errore quadratico medio Potenza media statistica dell’errore di misura: Nota bene:

7. Misura della matrice di covarianza Misura della matrice di covarianza del vettore degli osservati, utilizzando K realizzazioni del vettore y Coeff. di correlazione tra Yi e Yj Valore effettivo: 8444 4844 4484 4448 Misura_mat_cov = 8.0342 3.9898 4.1356 4.1333 3.9898 8.0060 3.9277 4.0460 4.1356 3.9277 8.1985 4.1425 4.1333 4.0460 4.1425 8.0706

8. File.m: segncost_matcov.m % Misura della matrice di covarianza degli osservati di segnale costante in rumore clcclear N=4; % numero di campioniK=5000; % numero di esperimenti (prove) Monte Carlo etas=10 % valor medio del segnale costantesig2s=4 % varianza del segnale costantesig2w=4 % varianza del rumore s=randn(1,K)*sqrt(sig2s)+etas; % genera ne valori di segnale costante con ddp Gaussianaw=randn(N,K)*sqrt(sig2w); % genera ne vettori di rumore con ddp Gaussiana, % organizzati in una matricey=repmat(s,N,1)+w; % calcola ne vettori osservati, organizzati in una matrice valvettmedio=mean(y.'); % calcola il vettor valor medio degli osservatiy=y-repmat(valvettmedio.',1,K); % sottrae dagli osservati il vettore valor mediovalmatcov=y*y.'/ne; % calcola la matrice di covarianza degli osservati,% sfruttando la loro organizzazione in matrice valmatcov % stampa la valutazione della matrice di covarianza % degli osservati

9. Simulazione stimatori e istogramma degli errori • Generare K=50000 realizzazioni del vettore osservato Y di N=2 campioni, per [istruzioni utili:randn, repmat] • Implementare i due stimatori, calcolare e visualizzare l’istogramma normalizzato dei loro errori di stima e confrontarli con gli andamenti teorici (ddp Gaussiane a valor medio nullo e varianza Pe ) [istruzioni utili: mean, hist, bar]

10. Istogrammi (normalizzati) degli errori K=50000 realizzazioni

11. File.m: stimesegncost_istog.m % Calcolo istogrammi degli errori degli algoritmi di misura del segnale costante clcclear N=2; % numero di campioniK=50000; % numero di esperimenti (prove) Monte Carlo etas=10 % valor medio del segnale costantesig2s=4 % varianza del segnale costantesig2w=4 % varianza del rumore g=sig2s/sig2w; % gamma s=randn(1,K)*sqrt(sig2s)+etas; % genera K valori di segnale costante con ddp Gaussianaw=randn(N,K)*sqrt(sig2w); % genera K vettori di rumore con ddp Gaussiana, % organizzati in una matricey=repmat(s,N,1)+w; % calcola K vettori osservati, organizzati in una matrice sa=mean(y); % calcola K misure di tipo asb=g*N/(g*N+1)*sa+etas/(g*N+1); % calcola K misure di tipo b esa=sa-s; % calcola i K errori delle misure di tipo aesb=sb-s; % calcola i K errori delle misure di tipo b [ista,aa]=hist(esa,50); % calcola l'istogramma degli errori delle misure di tipo a[istb,ab]=hist(esb,50); % calcola l'istogramma degli errori delle misure di tipo b figure % esegue i grafici degli istogrammi normalizzatibar(aa,ista/K/(aa(2)-aa(1)))xlabel('errore')ylabel('ddp')title('a')figurebar(ab,istb/K/(ab(2)-ab(1)))xlabel('errore')ylabel('ddp')title('b')

12. Curve di prestazioni degli algoritmi di misura • Ripetere la simulazione per un numero di campioni N variabile da 2 a 50, con 1000 realizzazioni, per , quindi valutare e visualizzare la potenza media degli errori Pe di entrambi gli algoritmi al variare di N [istruzioni utili:randn, repmat, for, mean, std, plot] • Ripetere la simulazione per variabile, per ; valutare e visualizzare la potenza media degli errori Pe di entrambi gli algoritmi al variare di g [istruzioni utili:randn, repmat, for, mean, std, plot]

13. K=1000 realizzazioni Prestazioni degli algoritmi in funzione di N • Il valore medio dell’errore è nullo per entrambi gli stimatori • La precisione della misura migliora con il numero N di campioni osservati • Il vantaggio dovuto alla conoscenza a priori si riduce all’aumentare di N

14. Al crescere di g per ed N fissati, l’utilità delle informazioni a priori si riduce Prestazioni degli algoritmi in funzione di g K=1000 realizzazioni • Il valore medio dell’errore è nullo per entrambi gli stimatori

15. File.m: stimesegncost_N.m % Calcolo della media e della varianza dell’errore dei due algoritmial variare del numero di campioni N clcclear Nmax=50; % numero massimo di campionine=1000; % numero di esperimenti (prove) Monte Carlo etas=10 % valor medio del segnale costantesig2s=4 % varianza del segnale costantesig2w=4 % varianza del rumore g=sig2s/sig2w; % gamma s=randn(1,ne)*sqrt(sig2s)+etas; % genera ne valori di segnale costante con d.d.p. Gaussianaw=randn(Nmax,ne)*sqrt(sig2w); % genera ne vettori di rumore con d.d.p. Gaussiana, % organizzati in una matricey=repmat(s,Nmax,1)+w; % calcola ne vettori osservati, organizzati in una matrice for indice=1:Nmax/2, % per numero di campioni N variabile a passi di 2: N=indice*2; sa=mean(y(1:N,:)); % calcola ne stime di tipo a con numero di campionisb=g*N/(g*N+1)*sa+etas/(g*N+1); % calcola ne stime di tipo b con numero di campioni N esa=sa-s; % calcola gli ne errori delle stime di tipo aesb=sb-s; % calcola gli ne errori delle stime di tipo bmsa(indice)=mean(sa); % calcola media delle stime di tipo avesa(indice)=mean(esa.^2); % calcola potenza errore delle stime di tipo amsb(indice)=mean(sb); % calcola media delle stime di tipo bvesb(indice)=mean(esb.^2); % calcola potenza errore delle stime di tipo b end a=2*[1:Nmax/2]; % esegue i graficifigureplot(a,msa,'r*',a,msb,'b*')xlabel('num. campioni')ylabel('medie')legend('stimatore a','stimatore b')axis([0 Nmax 0 2*etas])figureplot(a,vesa,'r*',a,vesb,'b*')hold on;plot(a,sig2w./a,'r',a,sig2w*g./(g*a+1),'b') % anche delle varianze teorichexlabel('num. campioni')ylabel('varianze')legend('stimatore a','stimatore b')

16. File.m: stimesegncost_gamma.m % Calcolo della media e della varianza dell’errore dei due algoritmial variare del rapporto gamma clcclear N=2; % numero di campionine=10000; % numero di esperimenti (prove) Monte Carlo etas=10 % valor medio del segnale costantesig2si=1 % varianza iniziale del segnale costantesig2w=4 % varianza del rumore for indice=1:20, % per rapporto gamma variabile: sig2s=sig2si*indice;g=sig2s/sig2w; % gamma s=randn(1,ne)*sqrt(sig2s)+etas; % genera ne valori di segnale costante con d.d.p. Gaussianaw=randn(N,ne)*sqrt(sig2w); % genera ne vettori di rumore con d.d.p. Gaussiana, % organizzati in una matricey=repmat(s,N,1)+w; % calcola ne vettori osservati, organizzati in una matrice sa=mean(y); % calcola ne stime di tipo asb=g*N/(g*N+1)*sa+etas/(g*N+1); % calcola ne stime di tipo b esa=sa-s; % calcola gli ne errori delle stime di tipo aesb=sb-s; % calcola gli ne errori delle stime di tipo bmsa(indice)=mean(sa); % calcola media delle stime di tipo avesa(indice)=mean(esa.^2);; % calcola potenza erroredelle stime di tipo amsb(indice)=mean(sb); % calcola media delle stime di tipo bvesb(indice)=mean(esb.^2);; % calcola potenza erroredelle stime di tipo b end a=sig2si*[1:20]/sig2w; % esegue i graficifigureplot(a,msa,'r*',a,msb,'b*')xlabel('gamma')ylabel('medie')legend('stimatore a','stimatore b')axis([min(a) max(a) 0 2*etas])figureplot(a,vesa,'r*',a,vesb,'b*')hold on;plot([a(1) a(20)],[sig2w/N sig2w/N],'r',a,sig2w*a./(a*N+1),'b') % anche delle varianze teoriche xlabel('gamma')ylabel('varianze')legend('stimatore a','stimatore b')title('sigw^2=4')