Part 2: Analysing the Spectral Content of an EEG Signal

This two-part series looks at the underlying mathematical theory and practical steps required to execute a spectral analysis of acquired electroencephalographic data. Spectral analysis, analysing how a signal’s power is distributed across a range of discrete frequencies, is fundamental to understanding an electroencephalogram as well as associated changes in psychological states. Solid, robust, replicable electrophysiology is built on the foundation of spectral analysis and one cannot emphasise enough the importance of the Fourier transform in this regard. The first part of this series deals with the general theory of the discrete Fourier transform (DFT), while the second part focuses on the mathematical/statistical properties of an EEG and how to apply a DFT in the context of biomedical signal-processing.

Main Image

Power spectral density of an electroencephalographic signal recorded over the occipital region during an eyes-open session. The shaded area around the spectrum is the 95% confidence interval around the mean of the estimate. The spectra has a number of prominent rhythmic components. The strong peak at 9Hz is the dominant alpha rhythm. One can also observe a strong beta component from 13Hz-17Hz. Together, these spectral features indicate that the person is alert, awake, and aware of their environment.

The software used for this two-part series is available for download here. The documentation for the software may be viewed here.
The Electroencephalogram, Random processes and the Power Spectral Density

EEG signals are seen as arising from random or stochastic processes, which means the observed/acquired signal, x[n], at any n is determined by a unique underlying probability distribution and is related across time by joint probability distributions (Oppenheim & Schafer, 2010, p. 1067-1068). For a random process, the power spectral density (PSD) is theoretically derived by taking the DFT of the autocorrelation sequence of the signal; however, in practical signal-processing, we can use the periodogram as an estimate of the PSD. Both methods are acceptable for purposes spectral estimation of random processes (ibid. Ch 10).

A signal generated by a random process is referred to as a realisation of that process. We may acquire multiple realisations of a given process during a single EEG recording. For analytical purposes, the acquired signal will be segmented into equal length segments—also referred to as epochs. For EEG signals, the segment length can be 0.5 – 2.5 seconds. Enforcing a segment length also ensures stationarity of the segment. Essentially, this means that the signal’s properties, likes its frequency content and amplitude, don’t change over time. In the context of random processes in general, and EEG signal-processing specifically, wide-sense stationarity is a requirement for application of the DFT to a signal or segment. Specifically, time-invariant spectral properties of the signal have a correspondence in relation to its statistical properties, such that:

  • the mean of the signal x[n] is time-invariant;
  • the autocovariance sequence is x[n] time-invariant, and
  • the variance of x[n] is finite across time points.

Enforcing a fixed length onto signal segments will also affect the lower bound of the analytical band-width, as at least one cycle of the lowest frequency would need to be represented in a segment to be amenable to analysis.

In EEG analysis, and in the analysis of most real-valued random signals, we often scale the power spectral estimate as the distribution of power in the frequency domain per unit frequency: Power/Hertz (\mu V^2/Hz), this scaled spectrum is specifically referred to as a power spectral density (PSD). In essence, it is a periodogram estimate; however, the scaling and normalisation applied differ from Eq 8. This is explained further in the following section.

Truncating and Scaling the PSD

Practically speaking, we are only interested in the values within the range from 0Hz – Nyquist frequency, where the Nyquist frequency is defined as half the sampling rate of the acquired data and is expressed in Hertz. As mentioned previously, the periodogram will provide spectral estimates up to and including the sampling rate of the data. Thus, if the sampling rate was 128Hz, we would only be interested in spectral estimates up to and including 64Hz. The reason is that any sampled data above the Nyquist frequency is subject to aliasing. Essentially, rhythms in the data occurring above the Nyquist frequency are under-sampled: one has not got enough datum-points to represent the entire rhythmic cycle of such components. So, if we were to attempt to reconstructed the original signal from the sampled data, the rhythms above the Nyquist would be reconstructed at an incorrect frequency, i.e. they appear as something they are not, an alias. So, we truncate V[k], removing all datum-points above the Nyquist frequency and denoted the new length as M:

    \[ 	V_T[k] = V[k], \quad 0 < k < M-1, \]


where the subscript T in V_T [k] denotes the truncated version of the DFT.

In order to conserve power correspondence between time and frequency domains—owing to the fact that spectral estimates above the Nyquist frequency were discarded, V_T [k] will be multiplied by 2 over the range 1<m<M-1, where m=0 is the DC-offset of the signal and m=M is the Nyquist frequency. For this particular estimate of the PSD, we compute a new window normalisation constant W:

    \[ 	W=LU=w[n]w[n]', \]


where ‘ denotes the vector transpose. In the case of a rectangular window W = 1 and can be substituted by window length L. In most practical instances, the scaled, normalised periodogram is taken as an estimate of the PSD:

    \[ 	I_{PSD}[k]=\left\{\begin{array}{@{}ll@{}} 		2 \triangle t\left(\frac{|V_T[k]|^2}{W}\right), & 1<k<M-1,\\ 		\;\, \triangle t\left(\frac{|V_T[k]|^2}{W}\right), & \text{otherwise}. 	\end{array}\right. \]


where the sampling interval \triangle t = 1/Fs; the frequency index of each bin, k, is given as, F= \frac{Fs}{N}k.

% compute power spectral density
W = w*w'; % normalisation constant
Vt = V(1:length(V)/2); % truncate DFT
It = zeros(1,length(Vt)); % zeropadded buffer
It(2:end-1) = 2*(1/Fs)*(abs(Vt(2:end-1)).^2/W); 
It(1) = 2*(1/Fs)*(abs(Vt(1)).^2/W);
It(end) = 2*(1/Fs)*(abs(Vt(end)).^2/W);
Perceptual Cycle
The power spectral density for a 100Hz sinusoidal embedded in white noise. The plot on the left shows the PSD. The plot on the right is the PSD in logarithmic scale.

Periodogram Averaging (Welch’s Periodogram)

Mentioned previously, for the purposes of spectral analysis using a DFT, x[n] will be segmented into K equal length, wide-sense stationary segments of length L—these maybe overlapping or not:

    \[ 	X_r[n] = x[rR+n]w[n], \quad 0 \leq n \geq L-1, \]


where x[n], 0 \le n \ge  Q-1 and Q is the length of the signal. The parameter, R, is a positive integer value. When R<  L, the segments will overlap; when R=L, the segments are non-overlapping (contiguous). The subscript r in x_r, is the index of the segment, 0 \ler \ge K-1. K is the total number of segments for which, (K-1)R+(L-1)\leQ-1. Thus, assuming R=L then K=Q/L (Oppenheim and Schafer, 2010, p. 868).

A windowed-DFT will be applied to each segment to yield an estimate, I_{PSD(r)}, computed as in Eq 11. The resultant local spectral estimates will be averaged to yield a final estimate:

    \[ 	\hat{I}_{PSD}[F]=\frac{1}{K}\sum_{r=0}^{k-1}I_{PSD(r)}[F],  \quad 0 \leq F \geq Nyquist. \]



% generate a signal (random noise)
rng default% reset random number generator

% signal
Fs = 512; % sampling rate
duration = 120; % 120 second duration signal
Q = Fs * duration;
x = 0 + 2.*randn(1,Q); % zero-mean random noise
t = 0:1/Fs:duration-1/Fs; % time-index (2 seconds)
x = cos(2*pi*100*t) + x; % sinusoidal at 100 Hz

% segmentation parameters
L = 1024;
R = 1024;
K = Q/L;
F = (Fs/L)*(0:L-1); % frequency index in Hertz

% compute periodogram
[w] = compHamming(L); % hamming window
PSDs = zeros(K,L/2);
 for r = 0:K-1
        % window the signal
        startInd = (r*R)+1;
        endInd = (startInd+L)-1;
        xr = x(1,startInd:endInd);
        v = xr.*w;
        [It] = computeScaledPSD(v,w,Fs);  % PSD
        PSDs(r+1,:) = It;

% average across segments
estPSD = mean(PSDs,1);
Perceptual Cycle
PSD’s of a 100Hz sinusoidal embedded in white noise (Q = 61440, R = 1024, L=1024, Fs = 512, Hamming window applied). The left panel displays a single realisation; the right panel displays a PSD’s of averaged across K = 60 realisations. Notice the large reduction in variance of the estimate in the second plot.
Perceptual Cycle
PSD’s of single-channel EEG data at the mid-line, occipital channel. The participant has his/her eyes open (Q = 24001, R=400, L=400, Fs = 200, Hamming window applied). The dataset was recorded at the Division of Biomedical Engineering University of Cape Town 2010. The plot on the left is a single realisation. The plot on the right is Welch’s periodogram estimate, K = 60. Note the prominent rhythmic component at 9Hz, this is the alpha rhythm. Also note the peak at 50Hz, which is the electrical noise from the power grid.
Perceptual Cycle
The top left panel shows the PSD (K = 240, Q = 24001, R = 100, L = 100, Fs = 200, Hamming window applied), here we see a large bias in the estimate which is most clear at 50Hz, where lines-noise from the power grid can be observed. An unbiased estimate would represent the component as a straight vertical line at 50Hz. On the plot we see the component much wider, indicating a large estimation bias. However, due to the short window length relative to total signal length, many segments are averaged, leading to an appreciable decrease in variance. The top right panel shows a PSD of the same signal (K = 1, L = Q = 24001), in this PSD we observe reduced bias, the 50Hz component is almost a straight vertical line at 50Hz; however, since there is only one segment, no averaging takes place and we observed a highly variable estimate. The segment would also not adhere to the requirement of wide-sense stationarity. The bottom panel displays a PSD of the same signal (K = 60, R = 400, L = 400), here we find an acceptable trade-off between estimation bias and variance. Also we ensure that the windowed signal is with in a range of what we would consider wide-sense stationary for a eyes-open EEG recording.

Baseline Standardisation

The PSD is represented in log or linear scale, or as a baseline standardised measure: measured relative to a baseline or reference period. The latter is often expressed as a percentage change in power, a z-score, or a difference score taken between logarithmically-scaled spectral estimates. The choice is dependent upon the application (Delorme & Makeig, 2004; Gerhold, 2017; Klimesch, 1999; Pfurtscheller & Lopes da Silva, 1999). Baseline standardised EEG measurements are used extensively in consumer neuroscience applications
(see Vecchiato, Cherubino, Trettel & Bablioni, 2013 for detailed methodologies applied to television advertising assessments using many of the above-mentioned baseline standardisation methodologies). The reader may also consult Gerhold (2017, Ch2, p. 11-54), wherein the author(s) discuss the various trade-offs involved in parametrisation of the DFT algorithm for analysis of event-related EEG data.

Measures of Band-Power

As Parceval’s theorem holds (which is generally stated for a continuous case of the Fourier transform), integrating values between any two points on a frequency domain function will yield energy in a sub-band bounded by two points. In the case of a PSD, summing values bounded by two frequency points f(1) and f(2) yields total band-power in that sub-band:

    \[ 	TotalBandPower=\sum_{F=f(1)}^{f(2)}  \hat{I}_{PSD}[F]. \]


Alternatively, the values including and between the upper and lower bound may be averaged to provide average band-power, a ubiquitous measure of oscillatory power of an EEG signal:

    \[ 	AvgBandPower=\frac{TotalBandPower}{J}, \]


where J is the total number of frequency bins summed across.

% compute band power measures
f1 = 7.5; f2 = 13;
f1 = find(F==f1); f2 = find(F==f2);
totalBandPower = sum(10*log10(estPSD(1,f1:f2)));
J = (f2-f1)+1; avgBandPower = totalBandPower/J;


Electricity Reflects the Content of the Mind

Through studying the power of the EEG signals in different sub-bands and how these band-power measures covary with task demands, in terms of cognitive or affective processes, we can gain an understanding of how the human brain (specifically how the cerebral cortex) responds to and handles information and even make predictions about populations of participants and consumers. For pivotal methodological and theoretical literature, see Klimesch (1999) and Klimesch, Sauseng, & Hanslmayr (2007). For consumer neuroscience applications see Ramsøy, Skov, Christensen, & Stahlhut (2018), Shestyuk, Kasinathan, Karapoondinott, Knight, & Gurumoorthy (2019) and Vecchiato, Cherubino, Trettel, & Babiloni (2013).

Delorme, A., & Makeig, S. (2004). EEGLAB: an open source toolbox for analysis of single-trial EEG dynamics including independent component analysis. Journal of Neuroscience Methods, 134, 9–21.

Gerhold, M. (2017). A study of event-related electrocortical oscillatory dynamics associated with cued motor-response inhibition during performance of the Go/NoGo task within a sample of prenatally alcohol-exposed children and age-matched controls. Cape Town: University of Cape Town.

Klimesch, W. (1999). EEG alpha and theta oscillations reflect cognitive and memory performance: a review and analysis. Brain Research Reviews, 29, 169–195.

Oppenheim, A. V., & Schafer, R. W. (2010). Discrete-Time Signal Processing (International Edition). Upper Saddle River: Pearson.

Pfurtscheller, G., & Lopes da Silva, F. H. (1999). Event-related EEG/MEG synchronization and desynchronization: basic principles. Clinical Neurophysiology, 110, 1842–1857.

Ramsøy, T. Z., Skov, M., Christensen, M. K., & Stahlhut, C. (2018). Frontal Brain Asymmetry and Willingness to Pay. Frontiers in Neuroscience, 12.

Shestyuk, A. Y., Kasinathan, K., Karapoondinott, V., Knight, R. T., & Gurumoorthy, R. (2019). Individual EEG measures of attention, memory, and motivation predict population level TV viewership and Twitter engagement. PLOS ONE, 14, e0214507.

Vecchiato, G., Cherubino, P., Trettel, A., & Babiloni, F. (2013). Neuroelectrical Brain Imaging Tools for the Study of the Efficacy of TV Advertising Stimuli and their Application to Neuromarketing. In Biosystems & Biorobotics. Retrieved from

This is an image

Published under: Consumer Neuroscience


Dr Matthew Gerhold is a cognitive-neuroscientist with a specialisation in behavioural neurosciences and non-invasive electrophysiological measurements. He has a wealth of experience extending over 14 years in different fields of the neurosciences: behavioural, cognitive, clinical and consumer. Gerhold is a doctoral graduate in the Health Sciences, conducting his research through the Division of Biomedical Engineering, Department of Human Biology at the University of Cape Town