Dr. Gari D. Clifford, University Lecturer & Associate Director,

1. Introduction to time-frequency analysisCentre for Doctoral Training in Healthcare Innovation Dr. Gari D. Clifford, University Lecturer & Associate Director, Centre for Doctoral Training in Healthcare Innovation, Institute of Biomedical Engineering, University of Oxford

2. Dealing with nonstationarities Physiological signals are dynamic; they exhibit time-varying statistics in both the time and frequency domain. weak-sense (or wide-sense stationarity or covariance stationarity) - mean and variance do not change significantly over a given window Strong/strictly stationary process - joint probability distribution does not change significantly when shifted in time or space. Higher order moments don’t change. Fourier analysis assumes weak stationarity What if the signal is not stationary? Physiological signals also exhibit activity that spans a range of time scales: E.g. The ECG signal contain three components with different time-scales: atrial depolarization represented by the P wave has a duration of 0.1 to 0.15 seconds; ventricular depolarization represented by the QRS complex 0.1 s; ventricular repolarization - T-wave 0.2 to 0.4 seconds.

3. 24 hour cardiovascular activity The inter-beat timing (average=heart rate), changes over much longer time scales; from minutes to hours, to days. Short term variations (1-10’s of s) - changes in the sympathetic (fight-and-flight) and parasympathetic (rest-and-digest) activity of the central nervous system acting on the heart. Longer term variations can be partially attributed to changes in activity and intrinsic circadian controls (that lead to sleep for example).

4. Tests for nonstationarity Parametric: Runs test (Bendat and Piersol 1986): Divide the series into time intervals of equal lengths. Compute a mean value (or other, see below) for each interval. Count the number of runs of mean values above and below the median value of the series. Compare the number of counts found to known probabilities of runs for random data and compare to look-up table. Non-parametric Sample auto-correlation function will neither cut off nor die down quickly, but rather will die down extremely slowly. (Bowerman and O'Connell 1979)

5. ‘Faking’ stationarity - Detrending Remove the mean or slope Remove a more complex shape – parabolic fitting Dangerous – but sometimes works! a(1)=rand(1,1); for(i=1:1000); a(i+1)=randn(1,1)+a(i); end b = detrend(a,'linear',[10:10:1000]);

6. ... but doesn’t always work Sum of sinusoids Up-Chirp n=1 (linear)

7. Where have we seen a chirp before?

8. Where have we seen a chirp before?

9. Why can’t we use the FFT?

10. Why can’t we use the FFT? Assumes stationarity! A sum of const frequency sinusoids So what can we do? Make shorter windows?

11. The STFT Short Term Fourier Transform

12. Spectrogram spectrogram(chirp,window,noverlap,[0:.1:1.2]',Fs) n=1 (linear) n=2 (quadratic)

13. STFT of VF • [S,F,T,P]=spectrogram(data,WINDOW,NOVERLAP,Fvect,Fs);

14. STFT of VF – longer time view [S,F,T,P]=spectrogram(data,WINDOW,NOVERLAP,Fvect,Fs);

15. VF – another STFT visualisation waterfall(P’)

16. STFT limitations – similar to FT A short analysis window results in good temporal resolution and poor spectral resolution. As the analysis window length increases its frequency resolution increases and its temporal resolution decreases.

17. STFT limitations The window w(t) determines the spectral and temporal resolution of the Short-Time Fourier Transform. Temporal resolution refers to the smallest time-separation below which two temporal events cannot be distinguished on the Spectrogram. Similarly, spectral resolution refers to the smallest frequency-separation below which two spectral events cannot be distinguished on the Spectrogram. A long window w(t) results in poor temporal resolution and good frequency resolution. Conversely, a short window w(t) results in good temporal resolution and poor spectral resolution. Time and frequency resolution can only be traded for one another, they cannot both be improved simultaneously.

18. Time-frequency trade-off The STFT limits our knowledge of signal activity to a two-dimensional cell in the time-frequency plane. The dimensions of the cell can be changed to favor either spectral or temporal resolution. The same temporal and spectral resolution applies to low and high frequency activity.

19. Revisit the Chirp Compare high temporal & high frequency resolution

20. Spike-train example Slow wave with periodic spike With STFT you can either resolve the location of the spikes or their frequency

21. How can we improve this trade-off? Is there a way to have good temporal and frequency resolution?

22. The Continuous Wavelet Transform A wavelet ψ (t)is shifted by a time τ and dilated or contracted by a factor a (the scale) prior to computing its correlation with the signal x(t). The correlation between the signal and the wavelet is defined as the integral of their product. The CWT maps x(t) into a bivariate function C(a,τ) that can be used to determine the similarity between x(t) and a wavelet scaled by a at given time. The correlation is localized in time, it is computed over an interval beginning at t = τ and ending t = τ + L where L is the duration of the wavelet. A time plot of the correlation C(a,τ) between the signal and the scaled wavelets is called a Scalogram. (Compare with spectrogram – scale is the inverse of frequency)

23. Construction of the scalogram Wavelet ψ (t)is shifted by a time τ and dilated or contracted by a factor a (the scale) prior to computing its correlation with the signal x(t). Coefficients of C(a,τ) represent the strength of the correlation at each scale at each point in time.

24. Time-Frequency Plane Tiling of CWT BUT - CWT automatically adjusts its time and frequency resolution depending on the scale of activity of interest by dilating or contracting the analysis window

25. CWT example Sinusoid With Frequency Breakdown STFT CWT

26. Chirp example with CWT CWT of Sinusoid With Linearly Increasing Frequency (n=1)

27. The Wavelet Transform The continuous wavelet transform (CWT) is a generalization of the STFT Allows for the analysis of non-stationary signals at multiple scales. Similar to the STFT, the CWT makes use of an analysis window to extract signal segments; Window is called a wavelet (or mother wavelet). Unlike the STFT, the analysis window or wavelet is not only convolved (translated), but dilated and contracted depending on the scale of activity under study. Wavelet dilation increases the CWT's sensitivity to long time-scale events Wavelet contraction increases its sensitivity to short time-scale events.

28. There are many types of wavelet basis functions In Matlab , type waveinfo('bior') etc

29. More wavelet basis functions

30. Properties of wavelets Similar to Fourier basis functions Dilations and translations of the Mother function or “analyzing wavelet” define our wavelet basis These basis functions are orthogonal – the correlation between different functions is zero Therefore we assume independence between the components They provide redundancy – and therefore error correction (1D →2D) They provide compact support – function is zero wherever it is not defined, vanishing at ∞ When zero distortion and anti-alias conditions are satisfied, perfect reconstruction is possible Basis functions are restricted in frequency – as we will now see when we consider the discrete wavelet transform Aside – we do not require windowing: why?

31. Discrete Wavelet Transform (DWT) The DWT can be viewed as a Filterbank STFT is an M-channel filterbank M parallel filters all with equal bandwidths but different center frequencies

32. Tree structured filterbank This filterbank architecture can be used to compute the coefficients of the continuous wavelet transform The coefficients are only computed over the dyadic scale (following a geometric sequence of ratio 2) Split into ‘details’ and ‘approximations’ at each stage H = filter transfer function a = approximations (not the scale here!) d = details

33. Details and approximations [lpdecon; hpdecon; lprecon; hprecon] = wfilters(‘bior4.4’); We can reconstruct original signal form all details and approximations

34. Wavelet decomposition ex 2 Four details and one approximation. Approximation is low freq trend

35. ECG filtering example ECG with mains noise ECG without mains noise With FIR filtering (& ringing) After wavelet filtering

36. Feature extraction from WT Wavelet ridges Wavelet modulus maxima

37. Applications Feature detection (ECG QT interval, spikes in EEG) Denoising Compression Trend extraction Further reading: Wavelet transforms and the ECG: a review Paul S Addison 2005 Physiol. Meas.26 R155