Ćwiczenia 4: Różnice pomiędzy wersjami

Z Brain-wiki
 
(Nie pokazano 15 pośrednich wersji utworzonych przez tego samego użytkownika)
Linia 3: Linia 3:
  
 
==Funkcja autokorelacji==
 
==Funkcja autokorelacji==
=== Korelacja i funkcja korelacji===
 
Dla przypomnienia: zagadnienia te poruszane były na [[Twierdzenia_o_splocie_i_o_próbkowaniu_(aliasing)#Korelacja_i_splot|wykładzie]]
 
 
Kowariancja:
 
 
::<math>\sigma _{xy} = \int {(x(t)-\bar{x}) (y(t)-\bar{y}) dt }</math>
 
 
Kowariancja dwóch zmiennych losowych jest miarą tego na ile dwie zmienne losowe mają podobne tendencje do zmian.
 
Przykład: w zmiennej <tt>x</tt> parami podajemy wartość pierwszej i drugiej zmiennej
 
<tt>
 
x = np.array([[0, 2], [1, 1], [2, 0]]).T
 
</tt>
 
W tej postaci łatwo zauważyć, że gdy pierwsza zmienna rośnie to druga maleje:
 
>>> x
 
array([[0, 1, 2],
 
        [2, 1, 0]])
 
Tę własność naszych zmiennych wyrażają elementy pozadiagonalne macierzy kowariancji:
 
>>> np.cov(x)
 
array([[ 1., -1.],
 
      [-1.,  1.]])
 
 
 
 
Współczynnik normalizacyjny:
 
 
::<math>\sigma _{ss} = \int {\left(s(t)-\bar{s}\right)^2 dt}</math>
 
Implementacja w Pythonie: [http://docs.scipy.org/doc/numpy/reference/generated/numpy.cov.html#numpy.cov numpy.cov]
 
 
 
Korelacja
 
::<math>\rho _{xy}= \frac{\sigma _{xy}}{\sqrt{\sigma _{xx} \sigma _{yy}}}</math>
 
 
Implementacja w Pythonie: [http://docs.scipy.org/doc/numpy/reference/generated/numpy.corrcoef.html numpy.corrcoef]
 
 
 
=== Funkcja korelacji wzajemnej===
 
=== Funkcja korelacji wzajemnej===
  
Dla sygnałów <math>x(t)</math> i <math>y(t)</math> sensowne jest rozszerzenie tych pojęć na wzajemne przesunięcia czasowe <math>\tau </math>:
 
  
::<math>\sigma _{xy}(\tau ) = \int {x(t) y(t+\tau ) dt} </math>
+
Dla sygnału  <math>x(t)</math>  możemy badać jak jest on do siebie podobny dla przesunięciaczasowego <math>\tau </math>:
a szczególnym przypadkiem jest funkcja autokorelacji.
 
Zwróćmy uwagę na związek funkcji korelacji ze splotem.
 
  
Implementacja w Pythonie: [http://docs.scipy.org/doc/numpy/reference/generated/numpy.correlate.html numpy.correlate]
+
::<math>\sigma _{xx}(\tau ) = \int {x(t) x(t+\tau ) dt} </math>
  
 +
Zwróćmy uwagę na związek funkcji korelacji:
 +
* z iloczynem skalarnym wektorów <math>x(t)</math> i jego przesuniętej wersji <math>x(t +\tau)</math>
 +
* ze splotem.
  
 
Dla skończonych dyskretnych sygnałów mamy estymatory korelacji:
 
Dla skończonych dyskretnych sygnałów mamy estymatory korelacji:
 +
::<math>R_{xx}(m) = E\lbrace x_{n+m}x^*_n \rbrace  = E\lbrace x_{n}x^*_{n-m} \rbrace </math>
  
::<math>R_{xy}(m) = E\lbrace x_{n+m}y^*_n \rbrace  = E\lbrace x_{n}y^*_{n-m} \rbrace </math>
+
oraz kowariancji:
 
+
::<math>C_{xx}(m) = E\lbrace (x_{n+m}-\bar{x}) (x_n-\bar{x})^*\rbrace  = R_{xx}(m) - \bar{x} \bar{x}^*</math>
kowariancji:
 
 
 
::<math>C_{xy}(m) = E\lbrace (x_{n+m}-\bar{x}) (y_n-\bar{y})^*\rbrace  = R_{xy}(m) - \bar{x} \bar{y}^*</math>
 
  
 
Funkcję <math>R</math> można estymować z jednej realizacji procesu (zakładamy jego ergodyczność):
 
Funkcję <math>R</math> można estymować z jednej realizacji procesu (zakładamy jego ergodyczność):
::<math> \widehat{R}_{xy}(m) = \left\lbrace  
+
::<math> \widehat{R}_{xx}(m) = \left\lbrace  
 
\begin{array}{ll}
 
\begin{array}{ll}
\sum _{n=0}^{N-m-1}{x_{n+m} y_n^*} & m \ge 0 \\
+
\sum _{n=0}^{N-m-1}{x_{n+m} x_n^*} & m \ge 0 \\
\widehat{R}_{yx}^*(-m) & m < 0
+
\widehat{R}_{xx}^*(-m) & m < 0
 
\end{array} \right.</math>
 
\end{array} \right.</math>
  
===Badanie zależności między sygnałami przy pomocy funkcji korelacji===  
+
===Zadanie: Jak obliczyć funkcję autokorelacji?===
 +
W pythonie mamy wydajną implementację funkcji autokorelacji (tak naprawdę to funkcji korelacji wzajemnej, której szczególnym przypadkiem jest autokorelacja)
 +
w: [http://docs.scipy.org/doc/numpy/reference/generated/numpy.correlate.html numpy.correlate]. Zanim jednak przejdziemy do rutynowego korzystania z tej funkcji warto sobie wyrobić pewną intuicję co do jej działania. Możemy sobie to działanie wyobrazić tak, że tworzymy dwie kopie naszego sygnału wydłużone zerami. jedna kopia ma oryginalny sygnał po środku, a w drugiej będzie on stopniowo przesuwany od pozycji takiej w której z nieruchomą kopią ma jeden punkt wspólny z lewej stron do takiej, kiedy ma jeden punkt wspólny z prawej strony. W każdym położeniu wzajemnym sygnałów musimy policzyć ich iloczyn skalarny.  Poniższy kod pomoże nam to zilustrować:
 +
 
 
<source lang = python>
 
<source lang = python>
def gabor(t0=0.5, sigma = 0.1, f = 10, T = 1, Fs = 128, phi =0 ):  
+
import numpy as np
dt = 1.0/Fs
+
import pylab as py
t = np.arange(0,T,dt)
+
 
s = np.exp( -0.5*((t-t0)/sigma)**2 )*np.cos(2*np.pi*f*t + phi)
+
# Niech nasz sygnał będzie:
return (s,t)
+
x = np.array([1,2,3,2,1])
</source>
+
 
# wygeneruj dwa sygnały długości <math>T=2s</math> próbkowane z częstością <math>f_s=128</math> Hz przy uzyciu funkcji <tt>gabor</tt>. Oba gabory mają częstość <math>f=10</math> Hz i <math>\sigma =0.1</math> s.  Oba sygnały <tt>s1</tt> i <tt>s2</tt> są centrowane na <math>t_0=0.5</math> s
+
# Ma on długość
# oblicz funkcję korelacji wzajemnej <tt>z = np.correlate(s1,s2,mode='full')</tt>
+
N = ...
# Jaka jest długość sygnału <tt>z</tt>?
 
# Wykreśl w funkcji odpowiednich skal czasu na dwóch subplotach: na górnym sygnały <tt>s1</tt> i <tt>s2</tt> a na dolnym <tt>z</tt>.  
 
# Zaobserwuj położenie maksimum funkcji korelacji wzajemnej. Jaki jest związek oscylacji w funkcji korelacji wzajemnej z oscylacjami funkcji <tt>s1</tt> i <tt>s2</tt>
 
  
Wskazówka: Związek między czasem <tt>t</tt> dla sygnałów <tt>s1</tt> i <tt>s2</tt> a skalą czasu dla korelacji <tt> f_corr_t</tt> można zapisać w Pythonie:
+
# przygotujmy miejsce na kopie sygnału
<source lang = python>
+
# najpierw kopia nieruchoma na środku
f_corr_t = np.zeros(2*len(t)-1)
+
x_2 = np.zeros(3*N-2, dtype = 'int')
f_corr_t[0:len(t)]= -t[len(t)::-1]
+
x_2[...] = x
f_corr_t[len(t):]=t[1:]
+
print('liczymy f. autokorelacji sygnału x= ',x)
 +
for k in range(-N,N-1):
 +
    print(' ')
 +
    x_1 = np.zeros(3*N-2, dtype = 'int')
 +
    x_1[...]=x  # to jest kopia sygnału x cofnięta o k próbek w stosunku do nieruchomej kopii
 +
    print('x1: \t'+str(x_1))
 +
    print('x2: \t'+str(x_2))
 +
    print(40*'-')
 +
    print('x1*x2\t', x_1*x_2)
 +
    f_corr[k+N] = ...
 +
    print('f_corr(',k+1,') = \t',...)
 +
# porównajmy wynik z wywołaniem funkcji bibliotecznej:
 +
...
 
</source>
 
</source>
  
+
===Zadanie: Funkcja auokorelacji funkcji okresowej===
*Powtórz punkty 1-5 zmieniając położenie sygnału <tt>s1</tt> od 0.5 do 0.1 z krokiem 0.1, oraz sygnału <tt>s2</tt> od 0.5 do 0.9 z krokiem 0.1.  Zaobserwuj związek między położeniem maksimum funkcji korelacji wzajemnej a odległością między centrami gaborów.
+
Funkcja okresowa jest do siebie podobna po przesunięciu o cały okres.  
 +
*Jak zatem wygląda funkcja autokorelacji fragmentu funkcji sinus?
 +
* Jak zmienia się ta funkcja wraz z długością fragmentu sygnału?
  
 
+
<source lang = python>
* Wykonaj analogiczne iteracje zachowując stałe położenie gaborów (dla obu <tt>t0 = 0.5</tt> zmieniaj natomiast częstość <tt>s2</tt> <math>f</math> od 10 Hz do 16 Hz co 2 Hz. Wymuś stały zakres osi y na -100:100(funkcja <tt>pylab.ylim((-100,100))</tt>)
 
*
 
<!--
 
 
import numpy as np
 
import numpy as np
 
import pylab as py
 
import pylab as py
  
def gabor(t0=0.5, sigma = 0.1, f = 10, T = 1, Fs = 128, phi =0 ):
+
Fs = 100
+
dt = 1/Fs
dt = 1.0/Fs
+
f = 10
t = np.arange(0,T,dt)
+
T = 0.2
s = np.exp( -0.5*((t-t0)/sigma)**2 )*np.cos(2*np.pi*f*t + phi)
+
t = np.arange(0,T,dt)
return (s,t)
+
s = np.sin(2*np.pi*f*t)
+
 
(s1,t) = gabor(t0=0.5,f=10,Fs=1000)
+
f_corr = ...# pełna wersja funkcji autokorelacji - za pomoca funkcji bibliotecznej
f_corr_t = np.zeros(2*len(t)-1)
+
tau  =...# wektor przesunięć wzajemnych - powinien mieć 0 na środku
f_corr_t[0:len(t)]= -t[len(t)::-1]
+
 
f_corr_t[len(t):]=t[1:]
+
py.subplot(2,1,1)
py.figure(1)
+
py.stem(t,s)
for i in range(4):
 
    (s1,t) = gabor(t0=0.5-i*0.1,f=10,Fs=1000) 
 
    (s2,t) = gabor(t0=0.5+i*0.1,f=10,Fs=1000)
 
    py.subplot(2,4,i+1)
 
    py.plot(t,s1,t,s2)
 
  
    py.subplot(2,4,4+(i+1))
+
py.subplot(2,1,2)
    z = np.correlate(s1,s2,mode='full')
+
py.stem(tau,f_corr)
    py.plot(f_corr_t,z)
 
py.figure(2)
 
for i in range(4):
 
    (s1,t) = gabor(t0=0.5,f=10,Fs=1000)
 
    f2 = 10+2*i
 
    (s2,t) = gabor(t0=0.5,f=f2,Fs=1000)
 
    py.subplot(2,4,i+1)
 
    py.plot(t,s1,t,s2)
 
    py.title(f2)
 
  
    py.subplot(2,4,4+(i+1))
 
    z = np.correlate(s1,s2,mode='full')
 
    py.plot(f_corr_t,z)
 
    py.ylim((-100,100))
 
 
py.show()
 
py.show()
-->
+
</source>
 +
===Zadanie: Funkcja auokorelacji białego szumuj===
 +
Analogicznie do zadania powyżej proszę zapoznać się z funknją autokorelacji dla białego szumu. Dla ustalenia uwagi niech to będą niezależne próbki z rozkładu N(0,1).
 +
 
 +
==Co musimy z tego zapamiętać ==
 +
* Jak liczymy funkcję autokorelacji?
 +
* Jak wygląda funkcja autokorelacji dla funkcji okresowej, a jak dla szumu?
 +
* To, że dla sygnałów deterministycznych możemy '''obliczyć''' funkcję autokorelacji, zaś dla procesów stochastycznych, które znamy jedynie na podstawie danych, musimy ją '''estymować'''.
 +
 
  
 
<!--
 
<!--

Aktualna wersja na dzień 11:13, 25 lis 2016

Analiza_sygnałów_-_ćwiczenia/AR_1


Funkcja autokorelacji

Funkcja korelacji wzajemnej

Dla sygnału [math]x(t)[/math] możemy badać jak jest on do siebie podobny dla przesunięciaczasowego [math]\tau [/math]:

[math]\sigma _{xx}(\tau ) = \int {x(t) x(t+\tau ) dt} [/math]

Zwróćmy uwagę na związek funkcji korelacji:

  • z iloczynem skalarnym wektorów [math]x(t)[/math] i jego przesuniętej wersji [math]x(t +\tau)[/math]
  • ze splotem.

Dla skończonych dyskretnych sygnałów mamy estymatory korelacji:

[math]R_{xx}(m) = E\lbrace x_{n+m}x^*_n \rbrace = E\lbrace x_{n}x^*_{n-m} \rbrace [/math]

oraz kowariancji:

[math]C_{xx}(m) = E\lbrace (x_{n+m}-\bar{x}) (x_n-\bar{x})^*\rbrace = R_{xx}(m) - \bar{x} \bar{x}^*[/math]

Funkcję [math]R[/math] można estymować z jednej realizacji procesu (zakładamy jego ergodyczność):

[math] \widehat{R}_{xx}(m) = \left\lbrace \begin{array}{ll} \sum _{n=0}^{N-m-1}{x_{n+m} x_n^*} & m \ge 0 \\ \widehat{R}_{xx}^*(-m) & m \lt 0 \end{array} \right.[/math]

Zadanie: Jak obliczyć funkcję autokorelacji?

W pythonie mamy wydajną implementację funkcji autokorelacji (tak naprawdę to funkcji korelacji wzajemnej, której szczególnym przypadkiem jest autokorelacja) w: numpy.correlate. Zanim jednak przejdziemy do rutynowego korzystania z tej funkcji warto sobie wyrobić pewną intuicję co do jej działania. Możemy sobie to działanie wyobrazić tak, że tworzymy dwie kopie naszego sygnału wydłużone zerami. jedna kopia ma oryginalny sygnał po środku, a w drugiej będzie on stopniowo przesuwany od pozycji takiej w której z nieruchomą kopią ma jeden punkt wspólny z lewej stron do takiej, kiedy ma jeden punkt wspólny z prawej strony. W każdym położeniu wzajemnym sygnałów musimy policzyć ich iloczyn skalarny. Poniższy kod pomoże nam to zilustrować:

import numpy as np
import pylab as py

# Niech nasz sygnał będzie:
x = np.array([1,2,3,2,1])

# Ma on długość
N = ...

# przygotujmy miejsce na kopie sygnału
# najpierw kopia nieruchoma na środku 
x_2 = np.zeros(3*N-2, dtype = 'int')
x_2[...] = x
print('liczymy f. autokorelacji sygnału x= ',x)
for k in range(-N,N-1):
    print(' ')
    x_1 = np.zeros(3*N-2, dtype = 'int')
    x_1[...]=x  # to jest kopia sygnału x cofnięta o k próbek w stosunku do nieruchomej kopii
    print('x1: \t'+str(x_1))
    print('x2: \t'+str(x_2))
    print(40*'-')
    print('x1*x2\t', x_1*x_2)
    f_corr[k+N] = ...
    print('f_corr(',k+1,') = \t',...)
# porównajmy wynik z wywołaniem funkcji bibliotecznej:
...

Zadanie: Funkcja auokorelacji funkcji okresowej

Funkcja okresowa jest do siebie podobna po przesunięciu o cały okres.

  • Jak zatem wygląda funkcja autokorelacji fragmentu funkcji sinus?
  • Jak zmienia się ta funkcja wraz z długością fragmentu sygnału?
import numpy as np
import pylab as py

Fs = 100
dt = 1/Fs
f = 10
T = 0.2
t = np.arange(0,T,dt)
s = np.sin(2*np.pi*f*t)

f_corr = ...# pełna wersja funkcji autokorelacji - za pomoca funkcji bibliotecznej
tau  =...# wektor przesunięć wzajemnych - powinien mieć 0 na środku

py.subplot(2,1,1)
py.stem(t,s)

py.subplot(2,1,2)
py.stem(tau,f_corr)

py.show()

Zadanie: Funkcja auokorelacji białego szumuj

Analogicznie do zadania powyżej proszę zapoznać się z funknją autokorelacji dla białego szumu. Dla ustalenia uwagi niech to będą niezależne próbki z rozkładu N(0,1).

Co musimy z tego zapamiętać

  • Jak liczymy funkcję autokorelacji?
  • Jak wygląda funkcja autokorelacji dla funkcji okresowej, a jak dla szumu?
  • To, że dla sygnałów deterministycznych możemy obliczyć funkcję autokorelacji, zaś dla procesów stochastycznych, które znamy jedynie na podstawie danych, musimy ją estymować.


Analiza_sygnałów_-_ćwiczenia/AR_1