% Example 3: Run GLMdenoise using finite impulse response (FIR) model
0001 %% Example 3: Run GLMdenoise using finite impulse response (FIR) model 0002 0003 %% Download dataset (if necessary) and add GLMdenoise to the MATLAB path 0004 0005 setup; 0006 0007 %% Load in the data 0008 0009 % Load in the data 0010 load('exampledataset.mat'); 0011 0012 % Check the workspace 0013 whos 0014 %% 0015 0016 %% Call GLMdenoisedata using FIR model 0017 0018 % It is typical to assume that the same fixed hemodynamic response function (HRF) 0019 % applies across all experimental conditions and all voxels. Indeed, this is the 0020 % default behavior of GLMdenoisedata. 0021 % 0022 % In contrast, the finite impulse response (FIR) model makes no assumption about the 0023 % shape of the hemodynamic response; it does this by estimating a separate parameter 0024 % for each time point in the response to each condition. 0025 % 0026 % Applying the FIR model can be useful for datasets where differences in timecourses 0027 % exist across conditions or voxels. Moreover, even if timecourse differences do not 0028 % necessarily exist, it may still be useful to apply the FIR model as a way of checking 0029 % that this is indeed the case. Also, applying the FIR model can be used as a means 0030 % of troubleshooting (e.g., it may reveal some error in the handling or preparation 0031 % of the data). 0032 % 0033 % GLMdenoisedata includes a flag that allows the FIR model to be used. Denoising 0034 % proceeds as usual because the FIR model simply changes the model of the experiment- 0035 % related effects in the data. 0036 0037 % Here we apply the FIR model to the example dataset. The 'fir' option indicates to 0038 % use the FIR model and the 20 indicates to estimate timecourses up to 20 time 0039 % points after condition onset. Notice that the design matrix passed to GLMdenoisedata 0040 % is in the usual format (a matrix that indicates condition onsets); there is no 0041 % need to construct the FIR basis functions as that is done internally by the code. 0042 % Also, notice that here we turn off bootstrapping to save computational time and 0043 % memory requirements. 0044 results = GLMdenoisedata(design,data,stimdur,tr, ... 0045 'fir',20,struct('numboots',0), ... 0046 'example3figures'); 0047 %% 0048 0049 % Let's inspect the results for an example voxel. By eye, it appears that the 0050 % FIR timecourses for this voxel are reasonably approximated by a single HRF. 0051 % Also, notice that the canonical HRF (red line) is substantially different from 0052 % the FIR timecourses, suggesting that there may be value in fitting an HRF to the data. 0053 xx = 53; yy = 27; zz = 1; 0054 tcs = squeeze(results.modelmd(xx,yy,zz,:,:)); % 35 conditions x 21 time points 0055 figure; hold on; 0056 plot(0:tr:20*tr,tcs'); % plot FIR timecourses 0057 plot(0:tr:20*tr,mean(tcs,1),'k-','LineWidth',3); % plot the mean timecourse 0058 hrf = getcanonicalhrf(stimdur,tr); % compute a canonical HRF 0059 % scale canonical HRF to best match the mean timecourse 0060 hrf = hrf * (pinv(hrf(1:21)')*mean(tcs,1)'); 0061 plot(0:tr:(length(hrf)-1)*tr,hrf,'r-','LineWidth',3); 0062 straightline(0,'h','k-'); 0063 ax = axis; 0064 axis([0 30 ax(3:4)]); 0065 xlabel('Time from condition onset (s)'); 0066 ylabel('BOLD signal (% change)'); 0067 title('FIR timecourses (black = mean, red = canonical HRF)'); 0068 %%