Merge remote-tracking branch 'Mark-Kramer/master'

Chapter2/Chapter 2 Solutions/Chapter_2_Solutions.txt → Chapter2/Chapter 2 Solutions/README.md

@@ -1,3 +1,3 @@
 Chapter 2 Solutions
 
-Each solution is a function that is either self-contained or makes a call to various helper functions within this folder. However, please ensure that Chapter 2’s files (particularly the .mat files) are within this folder or within your MATLAB working directory. If not, then the various commands load('Ch2-EEG-#.mat') will fail.
+Each solution is a function that is either self-contained or makes a call to various helper functions within this folder. However, please ensure that Chapter 2’s files (particularly the .mat files) are within this folder or within your MATLAB working directory. If not, then the various commands load('Ch2-EEG-#.mat') will fail.

Chapter3/Chapter 3 Solutions/Chapter_3_Analyze_EEG.m

@@ -0,0 +1,64 @@
+function Chapter_3_Analyze_EEG(EEG, t)
+ %INPUT:
+ % EEG: row vector of mean-normalized EEG data
+ % t: column vector of within trial timestamps
+
+ %Calculate Various Measures
+ sample_interval = t(2) - t(1);
+ sample_freq = 1/sample_interval;
+ T = max(t);
+ N = length(t);
+ freq_resolution = 1/T;
+ nyquist_freq = sample_freq/2;
+
+ %Plot Data
+ figure()
+ plot(t, EEG, 'k', 'LineWidth', 2)
+ xlabel('Time (seconds)')
+ ylabel('Potential (\muV)')
+ title('Original Data')
+ set(gca, 'FontSize', 14)
+
+ %Plot Biased Autocovariance Lags
+ [ac, lags] = xcorr(EEG, 100, 'biased');
+ figure()
+ plot(lags*sample_interval, ac, 'b', 'LineWidth', 2)
+ xlabel('Lag (seconds)')
+ ylabel('Autocovariance (biased)')
+ title('Autocovariance of the Data')
+ set(gca, 'FontSize', 14)
+
+ %Fourier Transform
+ xf = fft(EEG); %get fourier transform (complex double output)
+ Sxx = 2*sample_interval^2/T*(xf.*conj(xf)); %get spectrum
+ Sxx = Sxx(1:N/2+1); %ignore negative freqs b/c they are redundant (so you don't add them) 
+ Sxx(1) = 0; %set 0 frequency to 0 (usually ~ 10^-30 due to rounding errors) 
+ faxis = (0:freq_resolution:nyquist_freq);
+ figure()
+ plot(faxis, 10*log10(Sxx), 'r', 'LineWidth', 2)
+ xlim([0 120])
+ xlabel('Freq (Hz)')
+ ylabel('Power (dB)')
+ title('Periodogram of Data')
+ set(gca, 'FontSize', 14)
+
+ %Compute Spectrogram
+ w_size = round(sample_freq/2); %size of window of data you want to analyze
+ overlap = round(w_size*.95); %overlap when computing
+ nfft = w_size; % number of points to consider for each step
+ [S, F, T, P] = spectrogram(EEG, w_size, overlap, nfft, sample_freq);
+ % S: complex double fourier transform at row = freq col = time
+ % F: column vector of freqs
+ % T: row vector of times
+ % P: power ("\muV^2/Hz")
+ P(1, :) = 0; %set 0 frequency to 0 (usually ~ 10^-30 due to rounding errors) 
+ figure()
+ imagesc(T, F, 10*log10(P)) 
+ colorbar %add colorbar
+ axis xy %origin in lower left; otherwise imagesc flips it weird
+ ylim([0 70]) %constrain frequency to [0 70] Hz
+ xlabel('Time (seconds)')
+ ylabel('Freqs (Hz')
+ set(gca, 'FontSize', 16)
+
+end

Chapter3/Chapter 3 Solutions/Chapter_3_Question_1.m

@@ -0,0 +1,26 @@
+function Chapter_3_Question_1()
+ % Question 1 - Frequency Resolution
+
+ %Determine Frequency Resolution:
+ % df == sample_freq/num_data
+ % == sample_freq/(sample_freq * length_of_recording)
+ % == 1/length_of_recording
+ % therefore, length_of_recording = 1/df
+ desired_freq_resolution = 0.0001;
+ required_recording_length = 1/desired_freq_resolution;
+
+ %The necessary time for recording is 10,000 seconds (about 2hr45min). The
+ %disadvantages of this are:
+ % 1) Most experiments are stimulus-response experiments, so one stimulus
+ % every 2hr45 minutes is a bit too long.
+ % 2) Recording that much data would require massive amounts of storage
+ % and subsequent computation.
+ % 3) Electrodes drift over time, so the signal present in the first half
+ % hour of recording will be different than the signal present in the
+ % last half hour of recording (for example).
+ % 4) If you had 2hr45min of recording, it would be better to do a
+ % spectrogram analysis than a fine-tuned spectrum. In most
+ % experiments, a frequency resolution of 0.0001 Hz isn't useful
+ % anyway.
+
+end

Chapter3/Chapter 3 Solutions/Chapter_3_Question_2.m

@@ -0,0 +1,49 @@
+function Chapter_3_Question_2()
+ % Question 2 - Periodogram v. Spectrum
+
+ %Set Up
+ load('Ch3-EEG-1.mat')
+ dt = t(2) - t(1); %sample interval
+ Fs = 1/dt; %sample freq
+ df = 1/max(t); %frequency resolution
+ fNQ = Fs/2; %nyquist frequency
+ N = length(EEG);
+ T = max(t);
+ data = EEG - mean(EEG);
+
+ %FFT Itself (from the book chapter)
+ xf = fft(data); %get fourier transform (complex double output)
+ Sxx = 2*dt^2/T*(xf.*conj(xf)); %get spectrum
+ Sxx = Sxx(1:N/2+1); %ignore negative freqs b/c they are redundant (so you don't add them) 
+ faxis_Sxx = (0:df:fNQ); %create frequency axis
+
+ %Periodogram Inputs (to prevent defaults within function)
+ % window: a vector of equal length to EEG; explored in Chapter 4
+ % nfft: the number of points you would like to include in the analysis; explored in Chapter 4
+ % Fs: the sampling frequency
+ window = rectwin(N); %our above analysis used no window (or, essentially, a rectangular window) 
+ nfft = N; %our above analysis used all data points
+ [Pxx, faxis_Pxx] = periodogram(data, window, nfft, Fs);
+
+ %Set 0 Term To Zero
+ % it should be zero, but rounding errors in MATLAB cause output of ~ 10^-33 
+ Pxx(1) = 0;
+ Sxx(1) = 0;
+
+ %Plot Both to Show They Are The Same
+ figure()
+ subplot(2, 1, 1)
+ plot(faxis_Sxx, 10*log10(Sxx/max(Sxx)), 'k', 'LineWidth', 2)
+ xlabel('Frequency (Hz)')
+ ylabel('Power (dB)')
+ legend({'FFT Analysis (Chapter)'})
+ set(gca, 'FontSize', 14)
+
+ subplot(2, 1, 2)
+ plot(faxis_Pxx, 10*log10(Pxx/max(Sxx)), 'b', 'LineWidth', 2)
+ xlabel('Frequency (Hz)')
+ ylabel('Power (dB)')
+ legend({'FFT Analysis (Periodogram)'})
+ set(gca, 'FontSize', 14)
+
+end

Chapter3/Chapter 3 Solutions/Chapter_3_Question_3.m

@@ -0,0 +1,32 @@
+function Chapter_3_Question_3()
+ % Question 3 - Autocorrelation (Normalized Autocovariance)
+
+ %Set Up
+ load('Ch3-EEG-1.mat')
+ dt = t(2) - t(1);
+ data = EEG - mean(EEG);
+ var_data = var(data);
+
+ %Autocovariance (biased): 
+ [auto_cov, lags] = xcorr(data, 100, 'biased'); %biased autocovariance
+
+ %Autocorrelation (biased):
+ auto_corr = auto_cov/var_data;
+
+ %Plot
+ figure()
+ hold on
+ plot(lags.*dt, auto_cov, 'k', 'LineWidth', 2)
+ plot(lags.*dt, auto_corr, 'b', 'LineWidth', 2)
+ hold off
+ xlabel('Lag (seconds)')
+ ylabel('Autocorrelation or Autocovariance')
+ legend({'Autocovariance', 'Autocorrelation'})
+ set(gca, 'FontSize', 14)
+
+ % The autocovariance is an unnormalized measure of covariance. This can be
+ % normalized by taking the biased autocorrelations and dividing by the
+ % variance of the data. Overlaying the data (as I have above) shows this
+ % clear underlying relationship between both measures.
+
+end

Chapter3/Chapter 3 Solutions/Chapter_3_Question_4.m

@@ -0,0 +1,32 @@
+function Chapter_3_Question_4()
+ % Question 4 - Explore/Analyze Ch3-EEG-2.mat
+
+ %Set Up
+ load('Ch3-EEG-2.mat')
+ data = EEG - mean(EEG);
+
+ %Run Analysis (with helper function)
+ Chapter_3_Analyze_EEG(data, t)
+ % a) The sampling interval is: 0.001 seconds.
+ % The total duration of recording is 1 second.
+ % The frequency resolution is 1 Hz.
+ % The Nyquist frequency is 500 Hz.
+ % b) The data appears to be dominated by a clear single rhythmic
+ % component. The rhythm looks to be roughly 60 Hz.
+ % c) The autocovariance shows again a clear 60 Hz rhythm.
+ % d) The periodogram shows a hidden rhythm at roughly 30 Hz and 40 Hz.
+ % This was not visible from the EEG raw data. I am using the dB scale
+ % to more clearly show small fluctuations that might be obscured by
+ % the large 60 Hz oscillation. The choice to normalize by the max 
+ % power is merely one convention and, when plotting, obscures subsequent
+ % problems' time-dependent changes. Thus, I chose not to divide by the
+ % max power.
+ % e) The spectrogram shows that the rhythms do not change with time.
+ % f) The initial description of the data being dominated by a 60 Hz
+ % rhythm remains true. However, it does not show the simultaneous 30
+ % and 40 Hz rhythms that are ongoing. By computing the periodogram, we
+ % saw the presence of the 30 to 40 Hz rhythms. By computing the
+ % spectrogram, we saw that these rhythms remained stable throughout
+ % the entire period.
+
+end

Chapter3/Chapter 3 Solutions/Chapter_3_Question_5.m

@@ -0,0 +1,32 @@
+function Chapter_3_Question_5()
+ % Question 5 - Explore/Analyze Ch3-EEG-3.mat
+
+ %Set Up
+ load('Ch3-EEG-3.mat')
+ data = EEG - mean(EEG);
+
+ %Run Analysis (with helper function)
+ Chapter_3_Analyze_EEG(data, t)
+ % a) The sampling interval is: 0.002 seconds.
+ % The total duration of recording is 2 seconds.
+ % The frequency resolution is 0.5 Hz.
+ % The Nyquist frequency is 250 Hz.
+ % b) The data appears to be dominated by a 20 Hz rhythm at the beginning
+ % and end of the recording.The middle doesn't seem to have any
+ % rhythmic component.
+ % c) The autocovariance shows a clear 20 Hz rhythm.
+ % d) The periodogram shows the 20 Hz rhythm. The rest appears to just be
+ % noisy measurements of frequency, indicating no other clear rhythms
+ % present in the signal. I am using the dB scale to more clearly show
+ % small fluctuations. The choice to normalize by the max power is 
+ % merely one convention and, when plotting, obscures this problem's
+ % time-dependent changes. Thus, I chose not to divide by the max
+ % power.
+ % e) The spectrogram shows that the 20 Hz rhythm occurs at the beginning
+ % and end of the recording, but not during the middle portions. No
+ % other clear rhythms are seen.
+ % f) The initial description of the data being dominated by a 20 Hz
+ % rhythm at the beginning and end remains true. Subsequent analysis
+ % confirmed our original findings.
+
+end

Chapter3/Chapter 3 Solutions/Chapter_3_Question_6.m

@@ -0,0 +1,38 @@
+function Chapter_3_Question_6()
+ % Question 6 - Explore/Analyze Ch3-EEG-4.mat
+
+ %Set Up
+ load('Ch3-EEG-4.mat')
+ data = EEG - mean(EEG);
+
+ %Run Analysis (with helper function)
+ Chapter_3_Analyze_EEG(data, t)
+ % a) The sampling interval is: 0.0005 seconds.
+ % The total duration of recording is 4 seconds.
+ % The frequency resolution is 0.25 Hz.
+ % The Nyquist frequency is 1000 Hz.
+ % b) The data appears to start with a 15 Hz rhythm that, over time,
+ % transitions into a 35 or so Hz rhythm.
+ % c) The autocovariance shows a roughly 25 Hz rhythm. This actually makes
+ % sense, given that the data evolves over time and the autocorrelation
+ % assumes it doesn't. Thus, if my signal is going from 15 Hz to 35 Hz
+ % or so, then the average will be 25 Hz, the rhythm seen in the
+ % autocorrelation.
+ % d) The periodogram shows that all frequencies between 15 and 35 Hz are
+ % equally represented in the signal. The rest appears to just be
+ % noisy measurements of frequency, indicating no other clear rhythms
+ % present in the signal. I am using the dB scale to more clearly show
+ % small fluctuations. The choice to normalize by the max power is 
+ % merely one convention and, when plotting, obscures this problem's
+ % time-dependent changes. Thus, I chose not to divide by the max
+ % power.
+ % e) The spectrogram shows that the 15 Hz rhythm occurs at the beginning
+ % and, over time, develops into a 35 Hz rhythm. No other frequencies
+ % have obvious changes.
+ % f) The initial description of the data being changed from a 15 Hz
+ % rhythm to a 35 Hz rhythm holds true. However, the autocorrelation
+ % measure didn't quite pick this up, adding merit to the idea
+ % spectral analysis, specifically spectrograms, should be done when
+ % dealing with time-dependent changes in signal.
+
+end

Chapter3/Chapter 3 Solutions/Chapter_3_Question_7.m

@@ -0,0 +1,39 @@
+function Chapter_3_Question_7()
+ % Question 7 - Explore/Analyze Ch3-EEG-5.mat
+
+ %Set Up
+ load('Ch3-EEG-5.mat')
+ data = EEG - mean(EEG);
+
+ %Run Analysis (with helper function)
+ Chapter_3_Analyze_EEG(data, t)
+ % a) The sampling interval is: 0.001 seconds.
+ % The total duration of recording is 2 seconds.
+ % The frequency resolution is 0.5 Hz.
+ % The Nyquist frequency is 500 Hz.
+ % b) The data is odd, to say the least. There seems to be a somewhat 2 Hz
+ % square wave that a noisy mixture of other frequencies are riding on
+ % top of.
+ % c) The autocovariance shows a slope heading down as the lags increase.
+ % There are no apparent oscillations in the autocovariance plot.
+ % d) The periodogram shows that every 4 Hz there is a peak in oscillatory
+ % power. This is likely because a square wave can be fit by a series
+ % of related frequencies overlayed on top of one another.The rest 
+ % appears to just be noisy measurements of frequency, indicating no
+ % other clear rhythms present in the signal. I am using the dB scale
+ % to more clearly show small fluctuations. The choice to normalize 
+ % by the max power is merely one convention and, when plotting, 
+ % obscures this problem's time-dependent changes. Thus, I chose not
+ % to divide by the max power.
+ % e) The spectrogram shows a peak frequency roughly every 4 Hz at the
+ % beginning, then a brief pause (corresponding to the large vertical
+ % line in the data), then again a peak frequency roughly every 4 Hz.
+ % f) The initial description of the data has been confirmed by the
+ % spectrogram, but not really they autocovariance or the periodogram.
+ % This particular form of wave (square wave) is particularly difficult
+ % to interpret, due to all the overlying frequencies needed to
+ % generate it. While the spectrogram kind of shows it, the data itself
+ % is definitely more useful in terms of understanding what the signal
+ % is (a square wave with noise at the peaks and troughs).
+
+end

Chapter3/Chapter 3 Solutions/Chapter_3_Question_8.m

@@ -0,0 +1,32 @@
+function Chapter_3_Question_8()
+ % Question 8 - Sine Function Analysis
+
+ %Generate Data
+ dt = 0.001; %note typo in book: should be 0.001 second not 0.001 millisecond 
+ t = dt:dt:10;
+ data = sin(2 * pi .* (t.^2));
+
+ %Run Analysis (with helper function)
+ Chapter_3_Analyze_EEG(data, t)
+ % a) The sampling interval is: 0.001 seconds.
+ % The total duration of recording is 10 seconds.
+ % The frequency resolution is 0.1 Hz.
+ % The Nyquist frequency is 500 Hz.
+ % b) The data is much like Ch3-EEG-4.mat data in question 6. It begins at
+ % a low frequency and increases to a 20 Hz or so frequency by the end
+ % of the data.
+ % c) The autocovariance shows a 16 Hz or so rhythm.
+ % d) The periodogram shows an increase in the [0 20] frequency range. 
+ % I am using the dB scale to more clearly show small fluctuations. 
+ % The choice to normalize by the max power is merely one convention 
+ % and, when plotting, obscures this problem's time-dependent changes.
+ % Thus, I chose not to divide by the max power.
+ % e) The spectrogram shows a clear change in the oscillatory power that
+ % starts at 0 Hz and goes until about 20 Hz. This is exactly what we
+ % saw by looking at the data.
+ % f) The initial description of the data has been confirmed by the
+ % spectrogram, but not really they autocovariance or the periodogram.
+ % This is expected because our data was a sine wave whose frequency
+ % varied (increased) with time.
+
+end