-
Notifications
You must be signed in to change notification settings - Fork 9
Analysis using codes: example 2
Using iMap4, we created the single-trial 2D fixation duration map and smoothed at 1° of visual angle. Importantly, to keep in line with Miellet, et al (2012), spatial normalization was performed by Z-scoring the fixation map across all pixels independently for each trial (the result is identical without spatial normalization in this example). We also applied a mask generated with the default option.
%% Load data
clear all
clc
filename=dir('data*.mat');
% create condition table
sbj = repmat([1:30],1,4)';
group = repmat([ones(1,15) ones(1,15)*2],1,4)';
blinkspot= repmat(1:4,30,1);
blinkspot= blinkspot(:);
Tbl = dataset(sbj,group,blinkspot);
% deg0cau = 1:15;
% deg0as = 16:30;
% deg2cau = 31:45;
% deg2as = 46:60;
% deg5cau = 61:75;
% deg5as = 76:90;
% deg8cau = 91:105;
% deg8as = 106:120;
% parameters for smoothing
ySize = 382;
xSize = 390;
smoothingpic= 10;
[x, y] = meshgrid(-floor(xSize/2)+.5:floor(xSize/2)-.5, -floor(ySize/2)+.5:floor(ySize/2)-.5);
gaussienne = exp(- (x .^2 / smoothingpic ^2) - (y .^2 / smoothingpic ^2));
gaussienne = (gaussienne - min(gaussienne(:))) / (max(gaussienne(:)) - min(gaussienne(:)));
f_fil = fft2(gaussienne);
Nitem=length(filename);
Ntall=zeros(Nitem,1);
for item=1:Nitem
load(['data' num2str(item) '.mat'])
Trials=unique(summary(:,4));
Ntall(item)=length(Trials);
end
Trialall=sum(Ntall);
rawmapMatST = zeros(Trialall,ySize, xSize);
fixmapMatST = zeros(Trialall,ySize, xSize);
groupST=zeros(Trialall,1);
sbjST=zeros(Trialall,1);
blinkspotST=zeros(Trialall,1);
stDurST=zeros(Trialall,1);
ii=0;
for item=1:Nitem
load(['data' num2str(item) '.mat'])
Nfix=size(summary,1);
Trials=summary(:,4);
UTrials=unique(Trials);
Nt=length(UTrials);
% condition based
fixmaptmp=zeros(Nt,ySize, xSize);
rawmaptmp=zeros(Nt,ySize, xSize);
for it=1:Nt
ii=ii+1;
coordX=round(summary(Trials==UTrials(it),2));%switch matrix coordinate here
coordY=round(summary(Trials==UTrials(it),1));
intv=summary(Trials==UTrials(it),3);
indx1=coordX>0 & coordY>0 & coordX<ySize & coordY<xSize;
rawmap=full(sparse(coordX(indx1),coordY(indx1),intv(indx1),ySize,xSize));
f_mat = fft2(rawmap); % 2D fourrier transform on the points matrix
filtered_mat = f_mat .* f_fil;
smoothpic = real(fftshift(ifft2(filtered_mat)));
fixmaptmp(it,:,:)=(smoothpic-mean(smoothpic(:)))./std(smoothpic(:));%
fixmapMatST(ii,:,:)=fixmaptmp(it,:,:);
rawmaptmp(it,:,:)=rawmap;
rawmapMatST(ii,:,:)=rawmaptmp(it,:,:);
groupST(ii)=Tbl.group(item);
sbjST(ii)=Tbl.sbj(item);
blinkspotST(ii)=Tbl.blinkspot(item);
stDurST(ii)=sum(indx1);
end
if Nt<10
fixmapMat(item,:,:)=mean(fixmaptmp,1);
rawmapMat(item,:,:)=mean(rawmaptmp,1);
else
fixmapMat(item,:,:)=trimmean(fixmaptmp,30);
rawmapMat(item,:,:)=trimmean(rawmaptmp,30);
end
end
TblST=dataset(sbjST,groupST,blinkspotST,stDurST);
TblST.blinkspotST=nominal(TblST.blinkspotST);
TblST.groupST=nominal(TblST.groupST);
TblST.groupST(TblST.groupST=='1')='WC';
TblST.groupST(TblST.groupST=='2')='EA';
TblST.groupST=nominal(TblST.groupST);
%% rescale
scale=150/mean([xSize,ySize]);
[ySize2,xSize2]=size(imresize(ones(ySize,xSize),scale));
fixmapMat2ST=zeros(size(TblST,1),ySize2,xSize2);
for it=1:size(TblST,1)
fixmapMat2ST(it,:,:)=imresize(squeeze(fixmapMatST(it,:,:)),scale);
end
%% mean map
figure('NumberTitle','off','Name','Mean fixation bias');
race=unique(TblST.groupST);
condi=unique(TblST.blinkspotST);
ii=0;
for ig=1:length(race)
for ipp=1:length(condi)
indxtmp=TblST.groupST==race(ig) & TblST.blinkspotST==condi(ipp);
TempMap=squeeze(mean(fixmapMat2ST(indxtmp,:,:),1));
% subplot(2,4,(ig-1)*4+ipp)
subplot(4,2,ipp*2-(2-ig))
ii=ii+1;betanew(ii,:,:)=TempMap;
imagesc(TempMap);colorbar
axis square, axis off,
title(char(race(ig)))
end
end
%
figure('NumberTitle','off','Name','all fixation bias');
subplot(1,2,1)
imagesc(squeeze(mean(fixmapMat2ST,1)));axis('equal','off')
subplot(1,2,2)
masktmpST=squeeze(mean(fixmapMat2ST,1))>.0045;
imagesc(masktmpST);axis('equal','off')
We then applied a full model on the single-trial fixation duration map made used of the “single-trial” option in iMap4:
Pixel_Intensity ~ Observer culture + Blindspot size + Observer culture * Blindspot size + (1 | subject)
Only the predictor of subject was treated as random effects and the model was fitted with maximum likelihood estimation (ML).
%% imapLMM
tic
opt.singlepredi=1;
[LMMmap,lmexample]=imapLMM(fixmapMat2ST,TblST,masktmpST,opt, ...
'PixelIntensity ~ groupST + blinkspotST + groupST:blinkspotST + (1|sbjST)','DummyVarCoding','effect');
toc
First check the model fitting parameters:
%% plot model fitting
close all
opt1.type='model';
% perform contrast
[StatMap]=imapLMMcontrast(LMMmap,opt1);
% output figure;
imapLMMdisplay(StatMap,0)
Then we can performed ANOVA on the LMM
%% plot fixed effec(anova result)
% close all
opt=struct;% clear structure
mccopt=struct;
opt.type='fixed';
opt.alpha=.05;
mccopt.methods='bootstrap';
mccopt.bootgroup={'groupST'};
% mccopt.methods='FDR';
mccopt.nboot=1000;
% mccopt.permute=1;
mccopt.bootopt=1;
% mccopt.sigma=smoothingpic*scale;
mccopt.tfce=0;
% perform contrast
[StatMap]=imapLMMcontrast(LMMmap,opt);
% mulitple comparison correction
[StatMap_c]=imapLMMmcc(StatMap,LMMmap,mccopt,fixmapMat2ST);
% output figure;
imapLMMdisplay(StatMap_c,1)
We can performed other linear contrast:
%% Compute linear contrast (reproduce figure 2 as in the orignial paper)
% close all
opt=struct;% clear structure
mccopt=opt;
opt.type='predictor beta';
opt.alpha=.05;
opt.c={[-1 0 0 0 1 0 0 0];[0 -1 0 0 0 1 0 0];[0 0 -1 0 0 0 1 0];[0 0 0 -1 0 0 0 1];[1 0 0 -1 0 0 0 0];[0 0 0 0 1 0 0 -1]};
opt.name={'WC-EA NV';'WC-EA 2dg';'WC-EA 5dg';'WC-EA 8dg';'WC NV-8dg'; 'EA NV-8dg'};
% opt.c=limo_OrthogContrasts([3,2]);
% opt.name={'Spotlight';'Position';'Interaction'};
% opt.h={[0.005],[0.005],[0.005],[0.005],[0.005],[0.005]};
% opt.onetail='>';
% mccopt.methods='Randomfield';
mccopt.methods='bootstrap';
mccopt.bootgroup={'groupST'};
% mccopt.methods='FDR';
mccopt.nboot=1000;
% mccopt.permute=1;
mccopt.bootopt=1;
% mccopt.sigma=smoothingpic*scale;
mccopt.tfce=0;
% perform contrast
[StatMap]=imapLMMcontrast(LMMmap,opt);
[StatMap_c]=imapLMMmcc(StatMap,LMMmap,mccopt,fixmapMat2ST);
% output figure;
imapLMMdisplay(StatMap_c,1)
For each unique catigorical condition, we can compute an "above chance" fixation pattern map.
%% Single predictor (above chance fixate)
opt=struct;% clear structure
mccopt=opt;
opt.type='predictor beta';
opt.alpha=.05;
opt.c={[1 0 0 0 0 0 0 0];[0 1 0 0 0 0 0 0];[0 0 1 0 0 0 0 0];[0 0 0 1 0 0 0 0];[0 0 0 0 1 0 0 0];[0 0 0 0 0 1 0 0];[0 0 0 0 0 0 1 0];[0 0 0 0 0 0 0 1]};
opt.name={'EA-NV';'EA-2deg';'EA-5deg';'EA-8deg';'WC-NV';'WC-2deg';'WC-5deg';'WC-8deg'};
h0=mean(fixmapMat2ST(repmat(masktmpST,[size(fixmapMat2ST,1),1,1])==1));
opt.h={h0,h0,h0,h0,h0,h0,h0,h0};
% opt.c=limo_OrthogContrasts([3,2]);
% opt.name={'Spotlight';'Position';'Interaction'};
% opt.h={[0.005],[0.005],[0.005],[0.005],[0.005],[0.005]};
opt.onetail='>';
% mccopt.methods='Randomfield';
mccopt.methods='bootstrap';
mccopt.bootgroup={'groupST'};
% mccopt.methods='FDR';
mccopt.nboot=1000;
% mccopt.permute=1;
mccopt.bootopt=1;
% mccopt.sigma=smoothingpic*scale;
mccopt.tfce=0;
% perform contrast
[StatMap]=imapLMMcontrast(LMMmap,opt);
[StatMap_c]=imapLMMmcc(StatMap,LMMmap,mccopt,fixmapMat2ST);
% output figure;
imapLMMdisplay(StatMap_c,1,[],'parula')
This wiki is adapted from the original iMap4 guidebook.
If you have any questions about the iMap4 usage, please email [email protected]
Getting started
Theory
- Linear Mixed Models
- Pixel Wise Modeling and non-parametric statistics
- Family-wise error rate (FWER) under H0
- Power analysis of iMap4
Data structures and function usage
- Core functions
- Input Matrix
- LMMmap
- StatMap, Posthoc and figure outputs
- Other useful features and function
Example 1 (GUI)
- Background of Example 1
- Using the GUI (1): Import Data and label columns
- Using the GUI (2): Parameters and Conditions
- Using the GUI (3): Create smoothed fixation matrix
- Using the GUI (4): Optional for preprocessing
- Using the GUI (5): Descriptive Statistics Report
- Using the GUI (6): Spatial Mapping Using Linear Mixed Models
- Using the GUI (7): Hypothesis testing and Display results
- Using the GUI (8): Post-hoc analysis
Example 2 (Code)
Example 3 (Code)
Example 4 (Code)
Future development
Additional information