6
ECE version 6342 Fall 2008 MATLAB Project
Reminder: This assignment is an individual assignment, and you must NOT work closely with other students. You must design, write, test, and debug the program alone. Discussions with other students must be limited to understanding concepts taught in the class or understanding MATLAB, i.e., general stuff, not solutions specific to this project. Similarly, using any part of a project solution from a previous semester is not allowed. Failures to observe this policy will be pursued as academic honesty violations.
MEG Analysis
Goal
In this project you will analyze several files, each containing a set of MEG (magnetoencephalography) event-related field potentials. In doing so, you will practice numerous signal processing skills.
Program
Create an m-file named zxy.m, where z, x and y are your initials, last name first. For example, if your name is Alice Rogers Smith, the file would be sar.m. This file is your main program. In the steps below you are asked to create additional functions in .m files, similarly named with your initials: zxy_bpf_butter_lp.m, zxy_erp_average.m, etc. In this program you will insert many consecutive analysis steps, each creating a plot or other output. After each plot or other analysis output, execute a pause statement to allow the user to consider the results. Important: For ALL plots, the horizontal axis must be properly labeled, with correct unit values indicated (appropriate time or frequency values). When you plot a DFT, label the axis with units in Hertz, not radians. When you plot a time-domain signal, label the axis in milliseconds. The entire plot should have a title that clearly indicates what it is.
Data
Your data can be found on the course web site in a data file named dataN.mat, where the N for your data can be found in the file data_assignments.txt. 1. The format of the file is a MATLAB .mat file. 2. The file contains 71,200 sample values representing 200 evoked potential trials, 356 data points per trial. The first 356 points are the first trial, the next 356 points are the second trial, and so on. At a sampling rate of 254 Hz, that is a trial signal length of approximately 1.402 seconds. 3. The first 51 points in each trial are pre-stimulus data. The stimulus is presented simultaneous with the 52nd data point, which therefore corresponds to the time t=0.
version 6
Analysis Steps
The data should be processed in steps similar to the following. 1. Load your data file. A new data vector chdat_red will contain the MEG data values from the file. 2. Divide the data into 200 separate signal vectors containing the 200 trials. 3. Plot the 200 trials all on the same axis. Be sure the vertical axis is correctly located at the t=0 point. 4. Remove the pre-stimulus offset in each trial a. Compute the average value of the pre-stimulus data points, points 1-51. b. Subtract that average from the entire trial (all 356 points). c. Use these trials in the remaining steps. 5. Ensemble average zxy_erp_average.m a. Compute the ensemble average of these (unfiltered) trials. b. Plot the ensemble average. It will be noisy, you but should see the shape of the evoked response. 6. Lowpass filter each trial zxy_lpf_butter_lp.m a. Design a 4th-order Butterworth LPF with cutoff frequency at 16 Hz. b. Plot the frequency response magnitude and phase. c. Plot the impulse response. d. Apply the filter to each of the trial to create filtered trials. Use MATLABs filter function. 7. Ensemble average zxy_erp_average.m a. Compute the ensemble average of the filtered trials. b. Plot the ensemble average together with the plot of the average of the unfiltered trials. It will look less noisy, and there will be a delay (between the unfiltered average and the filtered average) due to the filter. c. When you do this, visually estimate the delay and write a message (just use the disp command) to the screen giving your estimate of the delay. (Yes, the unfiltered average is noisy, but do your best.) 8. DFT of ensemble average a. Compute the 512-pt DFT and the 1024-pt DFT of the ensemble average; plot the magnitude of both. b. Compute the dominant frequency, the frequency corresponding to the highest magnitude of the DFT; print the result. c. Compute the % power in the standard MEG bands: delta (0-4 Hz), theta (4-8 Hz), alpha (8-12 Hz), beta (above 12 Hz); print the results. 9. Linear phase filtering a. Design a Parks-McClellan linear phase FIR LPF zxy_lpf_pm.m. To choose the length of this filter, inspect the plot of the Butterworth filters impulse response above, and estimate the point where that impulse response has died down close to zero. b. Repeat ensemble average.
version 6 c. Plot the ensemble average together with the previous two ensemble averages from above. d. When you do this, visually estimate the delay caused by this filter, and write a message (just use the disp command) to the screen giving your estimate of that delay. It should be the same as half the length of this fi...