Autokorelacja

Relacja do stacjonarności w węższym sensie

Proces stacjonarny w węższym sensie (ściśle stacjonarny) charakteryzuje się niezmiennością całego łącznego rozkładu prawdopodobieństwa przy przesunięciu w czasie. Jeżeli proces ten posiada skończoną wartość oczekiwaną i skończoną wariancję, to automatycznie jest również stacjonarny w szerszym sensie; wówczas jego autokorelacja zależy wyłącznie od opóźnienia ${\displaystyle k}$.

Istnieją jednak procesy ściśle stacjonarne mające rozkłady o tzw. grubych ogonach (np. rozkład Cauchy'ego), dla których wartość oczekiwana i wariancja nie istnieją – w ich przypadku klasyczny współczynnik autokorelacji nie może być zdefiniowany.

Sygnały ciągłe

(1) Dla sygnału ciągłego ${\displaystyle f(t)}$ o skończonej energii funkcję autokorelacji definiuje się jako całkę splotową, będącą funkcją przesunięcia ${\displaystyle \tau }$:

${\displaystyle R_{ff}(\tau )\equiv f(\tau )*{\overline {f}}(-\tau )=\int _{-\infty }^{\infty }f(t+\tau ){\overline {f}}(t)\,dt}$

lub równoważnie

${\displaystyle R_{ff}(\tau )\equiv f(\tau )*{\overline {f}}(-\tau )=\int _{-\infty }^{\infty }f(t){\overline {f}}(t-\tau )\,dt}$

gdzie:

• ${\displaystyle \tau \in \mathbb {R} }$ jest przesunięciem czasowym (opóźnieniem),
• ${\displaystyle *}$ oznacza splot funkcji,
• ${\displaystyle {\overline {f}}(t)}$ to sprzężenie zespolone sygnału (dla sygnałów rzeczywistych ${\displaystyle {\overline {f}}(t)=f(t)}$).

(2) Jeśli analizowany sygnał jest okresowy o okresie ${\displaystyle T}$ (przy czym ma skończoną moc, ale nieskończoną energię), autokorelację dla przesunięcia czasowego ${\displaystyle \tau }$ oblicza się jako średnią po jednym okresie:

${\displaystyle R_{ff}(\tau )={\frac {1}{T}}\int _{0}^{T}f(t)\cdot {\overline {f}}(t-\tau )\,dt}$

Sygnały dyskretne

Dla sygnału dyskretnego ${\displaystyle \{x_{n}\}}$ złożonego z ${\displaystyle N}$ próbek odległych w czasie o ${\displaystyle dt}$, estymatę funkcji autokorelacji dla przesunięcia o ${\displaystyle k}$ próbek (czyli przesuniętych w czasie o ${\displaystyle \Delta t=k\cdot dt}$, ${\displaystyle k=0,\pm 1,\pm 2,\dots }$ ) określa wzór:

${\displaystyle R_{xx}(k)={\frac {1}{N-k}}\sum _{n=0}^{N-k-1}x_{n}\cdot {\overline {x}}_{n+k}}$

Próbkowanie sygnału

Sygnał ciągły zostaje poddany procesowi próbkowania (dyskretyzacji) z okresem próbkowania ${\displaystyle dt=0{,}01}$ s, co definiuje częstotliwość próbkowania jako:

${\displaystyle f_{s}={\frac {1}{dt}}=100{\text{ Hz}}}$

Zgodnie z twierdzeniem Nyquista-Shannona, maksymalna częstotliwość, jaką można zarejestrować bez aliasingu, wynosi:

${\displaystyle f_{\text{max}}={\frac {f_{s}}{2}}=50{\text{ Hz}}}$

Generowany sygnał dyskretny składa się z komponentu użytecznego (sinusoidy o częstotliwości ${\displaystyle f=1{\text{ Hz}}}$) oraz szumu losowego:

${\displaystyle y[i]=\sin(2\pi fdt\cdot i)+n[i]}$,

gdzie ${\displaystyle n[i]}$ reprezentuje dyskretny szum biały o rozkładzie normalnym, średniej zero i wariancji ${\displaystyle \sigma ^{2}=4}$. Ponieważ amplituda sinusoidy wynosi ${\displaystyle 1}$ (co oznacza moc składowej zmiennej równą ${\displaystyle 0{,}5}$), a moc szumu wynosi ${\displaystyle 4}$, stosunek mocy sygnału do szumu (SNR) jest bardzo niski i wynosi ${\displaystyle 1:8}$ (${\displaystyle -9{\text{ dB}}}$).

Funkcja autokorelacji

Autokorelacja sygnału dyskretnego mierzy stopień powtarzalności i podobieństwa sygnału do samego siebie przesuniętego w czasie o opóźnienie wynoszące ${\displaystyle r}$ próbek. Teoretyczna funkcja autokorelacji dla sygnału ${\displaystyle y[i]}$ jest definiowana wzorem:

${\displaystyle R_{y}[r]=E\{y[i]\cdot y[i+r]\}}$

Ponieważ sygnał wyjściowy jest sumą sinusoidy ${\displaystyle s[i]}$ i szumu ${\displaystyle n[i]}$

${\displaystyle y[i]=s[i]+n[i]}$,

to otrzymamy:

${\displaystyle R_{y}[r]=E{\big \{}{\big (}s[i]+n[i]{\big )}\cdot {\big (}s[i+r]+n[i+r]{\big )}{\big \}}}$

Dzięki liniowości operatora wartości oczekiwanej

${\displaystyle E\{A+B\}=E\{A\}+E\{B\}}$

możemy zapisać ${\displaystyle R_{y}[r]}$ jako sumę czterech oddzielnych wartości oczekiwanych (wcześniej mnożąc odpowiednie wyrażenia w nawiasach):

${\displaystyle R_{y}[r]=E\{s[i]\cdot s[i+r]\}+E\{s[i]\cdot n[i+r]\}+E\{n[i]\cdot s[i+r]\}+E\{n[i]\cdot n[i+r]\}}$

Każdy z tych czterech składników pełni inną rolę fizyczną i matematyczną:

(1) ${\displaystyle E\{s[i]\cdot s[i+r]\}={\tfrac {1}{2}}\cos(2\pi f\cdot dt\cdot r)}$

– autokorelacja sygnału użytecznego; dla fali sinusoidalnej autokorelacja zachowuje okresowość i przyjmuje postać funkcji cosinus

(2) ${\displaystyle E\{s[i]\cdot n[i+r]\}=0\quad }$

oraz

${\displaystyle \quad E\{n[i]\cdot s[i+r]\}=0}$

– wyrazy korelacji wzajemnej sygnału użytecznego i szumu zerują się, gdyż losowy szum jest całkowicie niezależny od deterministycznej sinusoidy,

(3) ${\displaystyle E\{n[i]\cdot n[i+r]\}={\begin{cases}\sigma ^{2}&{\text{dla }}r=0\\0&{\text{dla }}r\neq 0\end{cases}}\quad }$

– autokorelacja szumu białego jest różna od zera tylko wtedy, gdy mnożymy próbkę przez tę samą próbkę, czyli dla ${\displaystyle r=0}$, zaś dla przesunięcia ${\displaystyle r\neq 0}$ autokorelacja szumu zeruje się (jest tak dlatego, że szum biały nie wykazuje pamięci i jego kolejne próbki są całkowicie niezależne od siebie)

Ostatecznie mamy:

${\displaystyle R_{y}[r]={\frac {1}{2}}\cos(2\pi f\,dt\cdot r)+{\begin{cases}\sigma ^{2}&{\text{dla }}r=0\\0&{\text{dla }}r\neq 0\end{cases}}}$

Wynik ten oznacza, że funkcja autokorelacji dla każdego ${\displaystyle r>0}$ jest funkcją okresowa, o okresie równym okresowi składowej sinusoidalnej sygnału, pozwalając na wyodrębnienie jej z szumu; wpływ szumu teoretycznie spada do zera.

Widmo Gęstości Mocy (PSD)

Drugą metodą detekcji jest transformacja sygnału z dziedziny czasu do dziedziny częstotliwości za pomocą Dyskretnej Transformaty Fouriera (ang. Discrete Fourier Transformation – DFT), realizowanej algorytmem szybkiej transformaty Fouriera (ang. Fast Fourier Transformation – FFT):

${\displaystyle Y[k]=\sum _{i=0}^{N-1}y[i]\cdot e^{-j{\frac {2\pi }{N}}ki}}$

Widmo gęstości mocy (ang. Power Spectral Density – PSD) definiuje rozkład mocy sygnału w dziedzinie częstotliwości:

${\displaystyle S_{y}[k]={\frac {|Y[k]|^{2}}{N\cdot f_{s}}}}$

Szum biały charakteryzuje się płaskim widmem (jednostajna gęstość mocy w całym pasmie), natomiast sygnał sinusoidalny skupia swoją energię, tworząc wyraźny, ostry prążek (pik) dokładnie dla częstotliwości ${\displaystyle f=1{\text{ Hz}}}$.

Realizacja w kodzie programistycznym

Program implementuje powyższą teorię w języku Python przy użyciu bibliotek NumPy (obliczenia wektorowe) oraz Matplotlib (wizualizacja).

for r in range(m):
    kor_yn[r] = np.dot(yn[:N-r], yn[r:]) / (N-r)
    kor_y[r]  = np.dot(y[:N-r],  y[r:])  / (N-r)