%% Data Supplement %%
% Data from: Woolnough, Forseth, Rollo and Tandon (2019) Uncovering the Functional Anatomy of the Insula during Speech, eLife
%%
clear; close all

%% Figure 2 : Broadband Gamma
% time: Vector of time points in milliseconds relative to onset of audible speech
% bga: Activation in % broadband gamma over baseline, {regions x tasks}(patients x time)

load('Figure2_BGA.mat')
regions = {'PI' 'AI' 'STG' 'CS' 'IFG'};
tasks = {'Reading' 'Naming_{Mono}' 'Naming_{Multi}' 'Listening'};
cmp={'r' 'g' [0 179 255]./255 [255 199 100]./255 [179 0 255]./255};

figure('color','w')
for i = 1:length(tasks)
    subplot(1,4,i)
    hold on
        for j = 1:length(regions)
        patch([time fliplr(time)],[nanmean(bga{j,i},1)+(nanstd(bga{j,i},[],1)./sqrt(size(bga{j,i},1))) fliplr(nanmean(bga{j,i},1)-(nanstd(bga{j,i},[],1)./sqrt(size(bga{j,i},1))))],cmp{j},'EdgeColor','none','FaceAlpha',0.3)
        end
    
    for j = 1:length(regions)
        plot(time,nanmean(bga{j,i},1),'Color',cmp{j})
    end
    axis([-1000 1200 -25 300])
    if i==1
        xlabel('Time (ms)')
        ylabel('% BGA')
        legend(regions)
        legend boxoff
    end
end

%% Figure 3
% coords: [x y z] coordinates of electrodes in pial space
% bga: Activation in % broadband gamma over baseline, {tasks x hemi}(elecs x time windows)
% zscore: Significance of activation
% thresh: FDR corrected q<0.05 threshold

load('Figure3.mat')
tasks = {'Reading' 'Naming_{Mono}' 'Naming_{Multi}' 'Praxis' 'Listening'};
c=linspace(1,0);
cmp = [c' c' ones(100,1)];

figure('color','w')
for e = 1:5 %Cycle tasks
    for h = 1:2 %Cycle hemispheres
        for t = 1:2 %Cycle time windows
            ind = ((e-1)*4) + ((t-1)*2) + h;
            subplot(5,4,ind)
            sig = (zscore{e,h}(:,t)>thresh(e,t));   %Electrodes greater than FDR signficance threshold
            hold on
            scatter3(coords{e,h}(sig,1),coords{e,h}(sig,2),coords{e,h}(sig,3),max([bga{e,h}(sig,t) 4*ones(size(bga{e,h}(sig,t)))],[],2),zscore{e,h}(sig,t),'MarkerFaceColor','flat')
            scatter3(coords{e,h}(~sig,1),coords{e,h}(~sig,2),coords{e,h}(~sig,3),4,'k','MarkerFaceColor','flat')
            
            if h==1
                view(-90,0)
                ylim([-10 35])
                zlim([-42 28])
            else
                view(90,0)
                ylim([-10 35])
                zlim([-42 28])
            end
            axis equal
            
            if h==1 && t==1
                zlabel(tasks{e})
            end
        end
    end
end
colormap(cmp)
%% Figure 5
% time: Vector of time points in milliseconds relative to onset of audible speech
% bga: Activation in % broadband gamma over baseline, {elecs x patient}(trials x time)

load('Figure5.mat')
figure('color','w')
for i = 1:4
    subplot(2,2,i)
    hold on
            patch([time fliplr(time)],[nanmean(bga{i},1)+(nanstd(bga{i},[],1)./sqrt(size(bga{i},1))) fliplr(nanmean(bga{i},1)-(nanstd(bga{i},[],1)./sqrt(size(bga{i},1))))],'b','EdgeColor','none','FaceAlpha',0.3)
    plot(time,mean(bga{i},1))
    axis([-1000 1200 -50 150])
    if i==3
        xlabel('Time (ms)')
        ylabel('% BGA')
    end
end


%% Figure 6
% coords: [x y z] coordinates of electrodes in inflated pial space
% bga: Difference in % broadband gamma between Reading - Listening
% zscore: Significance of difference between Reading - Listening
% thresh: FDR corrected q<0.05 threshold

load('Figure6.mat')
c=linspace(0.929,0.494);
d=linspace(0.694,0.184);
e=linspace(0.125,0.556);
cmp = [c' d' e'];

figure('color','w')
sigp = (zscore>thresh);   %Electrodes greater than FDR signficance threshold - Positive
sign = (zscore<-thresh);    %Electrodes greater than FDR signficance threshold - Negative
nsig = (abs(zscore)<thresh);
hold on
scatter3(coords(sigp,1),coords(sigp,2),coords(sigp,3),max([bga(sigp) 4*ones(size(bga(sigp)))],[],2),zscore(sigp),'MarkerFaceColor','flat')
scatter3(coords(sign,1),coords(sign,2),coords(sign,3),max([-bga(sign) 4*ones(size(bga(sign)))],[],2),zscore(sign),'MarkerFaceColor','flat')
scatter3(coords(nsig,1),coords(nsig,2),coords(nsig,3),4,'k','MarkerFaceColor','flat')

view(-90,0)
axis equal
colormap(cmp)