
%%%%%%%%%%%
% Spike Heights
AllCellsSpikeHeightStimRFON=[];
AllCellsSpikeHeightStimRFOFF=[];
mod_idx_SpikeHeightStim=[];
AllCellsSpikeHeightShamRFON=[];
AllCellsSpikeHeightShamRFOFF=[];
mod_idx_SpikeHeightSham=[];

for sesi=1:length(allVars)
    s_numTransient=[];
    for celli=1:length(allVars{sesi}.s)
    s_numTransient(celli)=length(find([allVars{sesi}.s{celli}]> 0));
    end
    criteria=[];
    criteria=find(s_numTransient>50 & ~ismember(1:size(allVars{sesi}.dff,2),allVars{sesi}.badcells));

    if allVars{sesi}.isStim==1
        AllCellsSpikeHeightStimRFON=[AllCellsSpikeHeightStimRFON; allVars{sesi}.SpikeHeightTrialStim(criteria)];
        AllCellsSpikeHeightStimRFOFF=[AllCellsSpikeHeightStimRFOFF; allVars{sesi}.SpikeHeightTrialNoStim(criteria)];
        mod_idx_SpikeHeightStim=[mod_idx_SpikeHeightStim;allVars{sesi}.mod_idx_SpikeHeight(criteria)];
    else
        AllCellsSpikeHeightShamRFON=[AllCellsSpikeHeightShamRFON; allVars{sesi}.SpikeHeightTrialStim(criteria)];
        AllCellsSpikeHeightShamRFOFF=[AllCellsSpikeHeightShamRFOFF; allVars{sesi}.SpikeHeightTrialNoStim(criteria)];
        mod_idx_SpikeHeightSham=[mod_idx_SpikeHeightSham;allVars{sesi}.mod_idx_SpikeHeight(criteria)];
    end
end

% for Stim sessions
figure();hold on;
scatter(AllCellsSpikeHeightStimRFOFF,AllCellsSpikeHeightStimRFON,20);
scatter(AllCellsSpikeHeightStimRFOFF(find(All_p_shuffl_modidx_SpikeHeight_Stim>0 & All_p_shuffl_modidx_SpikeHeight_Stim<0.05)),...
    AllCellsSpikeHeightStimRFON(find(All_p_shuffl_modidx_SpikeHeight_Stim>0 & All_p_shuffl_modidx_SpikeHeight_Stim<0.05)),20,'r');
line([0 max(max(xlim),max(ylim))],[0 max(max(xlim),max(ylim))],'linestyle',':', 'color','k')
slope=AllCellsSpikeHeightStimRFOFF\AllCellsSpikeHeightStimRFON;
plot(AllCellsSpikeHeightStimRFOFF,slope*AllCellsSpikeHeightStimRFOFF,'linestyle',':', 'color','r')
xlabel('Sum of Spike Heights for RF-OFF (a.u.)')
ylabel('Sum of Spike Heights for RF-ON (a.u.)')
title({['total number of traces Stim pardigm=' num2str(length(find(All_p_shuffl_modidx_SpikeHeight_Stim>0)))];...
    ['total num of responding cells=' num2str(length(find(All_p_shuffl_modidx_SpikeHeight_Stim>0 & All_p_shuffl_modidx_SpikeHeight_Stim<0.05)))];...
    ['Lin Reg: y=' num2str(slope) '*x']});
set(gca,'fontweight','bold','fontsize',11,'linewidth',1)
%inset
axes('Position',[.2 .7 .2 .2]); hold on;
box on
scatter(AllCellsSpikeHeightStimRFOFF,AllCellsSpikeHeightStimRFON,10);
scatter(AllCellsSpikeHeightStimRFOFF(find(All_p_shuffl_modidx_SpikeHeight_Stim>0 & All_p_shuffl_modidx_SpikeHeight_Stim<0.05)),...
    AllCellsSpikeHeightStimRFON(find(All_p_shuffl_modidx_SpikeHeight_Stim>0 & All_p_shuffl_modidx_SpikeHeight_Stim<0.05)),10,'r');
set(gca,'xscale','log','yscale','log')
line([5e-1 5e2],[5e-1 5e2],'linestyle',':', 'color','k')
xlim([5e-1,5e2]);ylim([5e-1,5e2]);
set(gcf,'position',[  292    42    950   900])
set(gca,'fontweight','bold','fontsize',11,'linewidth',1)
saveas(gca,fullfile(cd,['SummaryAllCellsSpikeHeight_Stim.pdf']))
saveas(gca,fullfile(cd,['SummaryAllCellsSpikeHeight_Stim.svg']))
saveas(gca,fullfile(cd,['SummaryAllCellsSpikeHeight_Stim.fig']))
close all

% for Sham sessions
figure();hold on;
scatter(AllCellsSpikeHeightShamRFOFF,AllCellsSpikeHeightShamRFON,20);
scatter(AllCellsSpikeHeightShamRFOFF(find(All_p_shuffl_modidx_SpikeHeight_Sham>0 & All_p_shuffl_modidx_SpikeHeight_Sham<0.05)),...
    AllCellsSpikeHeightShamRFON(find(All_p_shuffl_modidx_SpikeHeight_Sham>0 & All_p_shuffl_modidx_SpikeHeight_Sham<0.05)),20,'r');
line([0 max(max(xlim),max(ylim))],[0 max(max(xlim),max(ylim))],'linestyle',':', 'color','k')
slope=AllCellsSpikeHeightShamRFOFF\AllCellsSpikeHeightShamRFON;
plot(AllCellsSpikeHeightShamRFOFF,slope*AllCellsSpikeHeightShamRFOFF,'linestyle',':', 'color','r')
xlabel('Sum of Spike Heights for RF-OFF (a.u.)')
ylabel('Sum of Spike Heights for RF-ON (a.u.)')
title({['total number of traces Sham pardigm=' num2str(length(find(All_p_shuffl_modidx_SpikeHeight_Sham>0)))];...
    ['total num of responding cells=' num2str(length(find(All_p_shuffl_modidx_SpikeHeight_Sham>0 & All_p_shuffl_modidx_SpikeHeight_Sham<0.05)))];...
        ['Lin Reg: y=' num2str(slope) '*x']});
set(gca,'fontweight','bold','fontsize',11,'linewidth',1)
%inset
axes('Position',[.2 .7 .2 .2]); hold on;
box on
scatter(AllCellsSpikeHeightShamRFOFF,AllCellsSpikeHeightShamRFON,10);
scatter(AllCellsSpikeHeightShamRFOFF(find(All_p_shuffl_modidx_SpikeHeight_Sham>0 & All_p_shuffl_modidx_SpikeHeight_Sham<0.05)),...
    AllCellsSpikeHeightShamRFON(find(All_p_shuffl_modidx_SpikeHeight_Sham>0 & All_p_shuffl_modidx_SpikeHeight_Sham<0.05)),10,'r');
set(gca,'xscale','log','yscale','log')
line([5e-1 5e2],[5e-1 5e2],'linestyle',':', 'color','k')
xlim([5e-1,5e2]);ylim([5e-1,5e2]);
set(gcf,'position',[  292    42    950   900])
set(gca,'fontweight','bold','fontsize',11,'linewidth',1)
saveas(gca,fullfile(cd,['SummaryAllCellsSpikeHeight_Sham.pdf']))
saveas(gca,fullfile(cd,['SummaryAllCellsSpikeHeight_Sham.svg']))
saveas(gca,fullfile(cd,['SummaryAllCellsSpikeHeight_Sham.fig']))
close all

%compare distribution of mod-idx
p=ranksum(mod_idx_SpikeHeightStim(find(All_p_shuffl_modidx_SpikeHeight_Stim>0)),mod_idx_SpikeHeightSham(find(All_p_shuffl_modidx_SpikeHeight_Sham>0)));
fig=figure(); hold on;
histogram(mod_idx_SpikeHeightStim(find(All_p_shuffl_modidx_SpikeHeight_Stim>0)),60,'DisplayStyle','stairs','Normalization','probability','EdgeColor','r')
histogram(mod_idx_SpikeHeightSham(find(All_p_shuffl_modidx_SpikeHeight_Sham>0)),60,'DisplayStyle','stairs','Normalization','probability','EdgeColor','b')
line(0*[1 1],ylim,'linestyle',':','color','k')
xlim([-1 1])
xlabel('Modulation Index (a.u.)')
ylabel('count #')
legend({'Stim','Sham'},'box','off')
title({['SpikeHeight-Stim-vs-Sham'];['ranksum p=' num2str(p)]})
set(gca,'fontweight','bold','fontsize',11,'linewidth',1)
set(fig,'Renderer', 'painters')
saveas(gca,fullfile(cd,['CompareDistModIdx_SpikeHeightStimVsSham.pdf']))
saveas(gca,fullfile(cd,['CompareDistModIdx_SpikeHeightStimVsSham.svg']))
saveas(gca,fullfile(cd,['CompareDistModIdx_SpikeHeightStimVsSham.fig']))
close all

%%%%%%%%%%%
% Spike Events
AllCellsSpikeEventsStimRFON=[];
AllCellsSpikeEventsStimRFOFF=[];
mod_idx_SpikeEventsStim=[];
AllCellsSpikeEventsShamRFON=[];
AllCellsSpikeEventsShamRFOFF=[];
mod_idx_SpikeEventsSham=[];

for sesi=1:length(allVars)
    s_numTransient=[];
    for celli=1:length(allVars{sesi}.s)
    s_numTransient(celli)=length(find([allVars{sesi}.s{celli}]> 0));
    end
    criteria=[];
    criteria=find(s_numTransient>50 & ~ismember(1:size(allVars{sesi}.dff,2),allVars{sesi}.badcells));

    if allVars{sesi}.isStim==1
        AllCellsSpikeEventsStimRFON=[AllCellsSpikeEventsStimRFON; allVars{sesi}.NumEventsTrialStim(criteria)];
        AllCellsSpikeEventsStimRFOFF=[AllCellsSpikeEventsStimRFOFF; allVars{sesi}.NumEventsTrialNoStim(criteria)];
        mod_idx_SpikeEventsStim=[mod_idx_SpikeEventsStim;allVars{sesi}.mod_idx_NumEvents(criteria)];
    else
        AllCellsSpikeEventsShamRFON=[AllCellsSpikeEventsShamRFON; allVars{sesi}.NumEventsTrialStim(criteria)];
        AllCellsSpikeEventsShamRFOFF=[AllCellsSpikeEventsShamRFOFF; allVars{sesi}.NumEventsTrialNoStim(criteria)];
        mod_idx_SpikeEventsSham=[mod_idx_SpikeEventsSham;allVars{sesi}.mod_idx_NumEvents(criteria)];
    end
end

% for Stim sessions
figure();hold on;
scatter(AllCellsSpikeEventsStimRFOFF,AllCellsSpikeEventsStimRFON,20);
scatter(AllCellsSpikeEventsStimRFOFF(find(All_p_shuffl_modidx_SpikeEvnt_Stim>0 & All_p_shuffl_modidx_SpikeEvnt_Stim<0.05)),...
    AllCellsSpikeEventsStimRFON(find(All_p_shuffl_modidx_SpikeEvnt_Stim>0 & All_p_shuffl_modidx_SpikeEvnt_Stim<0.05)),20,'r');
line([0 max(max(xlim),max(ylim))],[0 max(max(xlim),max(ylim))],'linestyle',':', 'color','k')
slope=AllCellsSpikeEventsStimRFOFF\AllCellsSpikeEventsStimRFON;
plot(AllCellsSpikeEventsStimRFOFF,slope*AllCellsSpikeEventsStimRFOFF,'linestyle',':', 'color','r')
xlabel('Sum of Spike Events for RF-OFF (a.u.)')
ylabel('Sum of Spike Events for RF-ON (a.u.)')
title({['total number of traces Stim pardigm=' num2str(length(find(All_p_shuffl_modidx_SpikeEvnt_Stim>0)))];...
    ['total num of responding cells=' num2str(length(find(All_p_shuffl_modidx_SpikeEvnt_Stim>0 & All_p_shuffl_modidx_SpikeEvnt_Stim<0.05)))];...
    ['Lin Reg: y=' num2str(slope) '*x']});
set(gca,'fontweight','bold','fontsize',11,'linewidth',1)
%inset
axes('Position',[.2 .7 .2 .2]); hold on;
box on
scatter(AllCellsSpikeEventsStimRFOFF,AllCellsSpikeEventsStimRFON,10);
scatter(AllCellsSpikeEventsStimRFOFF(find(All_p_shuffl_modidx_SpikeEvnt_Stim>0 & All_p_shuffl_modidx_SpikeEvnt_Stim<0.05)),...
    AllCellsSpikeEventsStimRFON(find(All_p_shuffl_modidx_SpikeEvnt_Stim>0 & All_p_shuffl_modidx_SpikeEvnt_Stim<0.05)),10,'r');
set(gca,'xscale','log','yscale','log')
line([5e-1 5e2],[5e-1 5e2],'linestyle',':', 'color','k')
xlim([5e-1,5e2]);ylim([5e-1,5e2]);
set(gcf,'position',[  292    42    950   900])
set(gca,'fontweight','bold','fontsize',11,'linewidth',1)
saveas(gca,fullfile(cd,['SummaryAllCellsSpikeEvents_Stim.pdf']))
saveas(gca,fullfile(cd,['SummaryAllCellsSpikeEvents_Stim.svg']))
saveas(gca,fullfile(cd,['SummaryAllCellsSpikeEvents_Stim.fig']))
close all

% for Sham sessions
figure();hold on;
scatter(AllCellsSpikeEventsShamRFOFF,AllCellsSpikeEventsShamRFON,20);
scatter(AllCellsSpikeEventsShamRFOFF(find(All_p_shuffl_modidx_SpikeEvnt_Sham>0 & All_p_shuffl_modidx_SpikeEvnt_Sham<0.05)),...
    AllCellsSpikeEventsShamRFON(find(All_p_shuffl_modidx_SpikeEvnt_Sham>0 & All_p_shuffl_modidx_SpikeEvnt_Sham<0.05)),20,'r');
line([0 max(max(xlim),max(ylim))],[0 max(max(xlim),max(ylim))],'linestyle',':', 'color','k')
slope=AllCellsSpikeEventsShamRFOFF\AllCellsSpikeEventsShamRFON;
plot(AllCellsSpikeEventsShamRFOFF,slope*AllCellsSpikeEventsShamRFOFF,'linestyle',':', 'color','r')
xlabel('Sum of Spike Events for RF-OFF (a.u.)')
ylabel('Sum of Spike Events for RF-ON (a.u.)')
title({['total number of traces Sham pardigm=' num2str(length(find(All_p_shuffl_modidx_SpikeEvnt_Sham>0)))];...
    ['total num of responding cells=' num2str(length(find(All_p_shuffl_modidx_SpikeEvnt_Sham>0 & All_p_shuffl_modidx_SpikeEvnt_Sham<0.05)))];...
        ['Lin Reg: y=' num2str(slope) '*x']});
set(gca,'fontweight','bold','fontsize',11,'linewidth',1)
%inset
axes('Position',[.2 .7 .2 .2]); hold on;
box on
scatter(AllCellsSpikeEventsShamRFOFF,AllCellsSpikeEventsShamRFON,10);
scatter(AllCellsSpikeEventsShamRFOFF(find(All_p_shuffl_modidx_SpikeEvnt_Sham>0 & All_p_shuffl_modidx_SpikeEvnt_Sham<0.05)),...
    AllCellsSpikeEventsShamRFON(find(All_p_shuffl_modidx_SpikeEvnt_Sham>0 & All_p_shuffl_modidx_SpikeEvnt_Sham<0.05)),10,'r');
set(gca,'xscale','log','yscale','log')
line([5e0 5e3],[5e0 5e3],'linestyle',':', 'color','k')
xlim([5e0,5e3]);ylim([5e0,5e3]);
set(gcf,'position',[  292    42    950   900])
set(gca,'fontweight','bold','fontsize',11,'linewidth',1)
saveas(gca,fullfile(cd,['SummaryAllCellsSpikeEvents_Sham.pdf']))
saveas(gca,fullfile(cd,['SummaryAllCellsSpikeEvents_Sham.svg']))
saveas(gca,fullfile(cd,['SummaryAllCellsSpikeEvents_Sham.fig']))
close all

%compare distribution of mod-idx
p=ranksum(mod_idx_SpikeEventsStim(find(All_p_shuffl_modidx_SpikeEvnt_Stim>0)),mod_idx_SpikeEventsSham(find(All_p_shuffl_modidx_SpikeEvnt_Sham>0)));
fig=figure(); hold on;
histogram(mod_idx_SpikeEventsStim(find(All_p_shuffl_modidx_SpikeEvnt_Stim>0)),60,'DisplayStyle','stairs','Normalization','probability','EdgeColor','r')
histogram(mod_idx_SpikeEventsSham(find(All_p_shuffl_modidx_SpikeEvnt_Sham>0)),60,'DisplayStyle','stairs','Normalization','probability','EdgeColor','b')
line(0*[1 1],ylim,'linestyle',':','color','k')
xlim([-1 1])
xlabel('Modulation Index (a.u.)')
ylabel('count #')
legend({'Stim','Sham'},'box','off')
title({['SpikeEvents-Stim-vs-Sham'];['ranksum p=' num2str(p)]})
set(gca,'fontweight','bold','fontsize',11,'linewidth',1)
set(fig,'Renderer', 'painters')
saveas(gca,fullfile(cd,['CompareDistModIdx_SpikeEventsStimVsSham.pdf']))
saveas(gca,fullfile(cd,['CompareDistModIdx_SpikeEventsStimVsSham.svg']))
saveas(gca,fullfile(cd,['CompareDistModIdx_SpikeEventsStimVsSham.fig']))
close all



%%%%%%%%%%%
% Area under curve

mod_idx_AreaCurveStim=[];
mod_idx_AreaCurveSham=[];

for sesi=1:length(allVars)
    s_numTransient=[];
    for celli=1:length(allVars{sesi}.s)
    s_numTransient(celli)=length(find([allVars{sesi}.s{celli}]> 0));
    end
    criteria=[];
    criteria=find(s_numTransient>50 & ~ismember(1:size(allVars{sesi}.dff,2),allVars{sesi}.badcells));

    if allVars{sesi}.isStim==1
        mod_idx_AreaCurveStim=[mod_idx_AreaCurveStim;allVars{sesi}.mod_idx_areacurve(criteria)];
    else
        mod_idx_AreaCurveSham=[mod_idx_AreaCurveSham;allVars{sesi}.mod_idx_areacurve(criteria)];
    end
end

%compare distribution of mod-idx
p=ranksum(mod_idx_AreaCurveStim(find(All_p_shuffl_modidx_areacurv_Stim>0)),mod_idx_AreaCurveSham(find(All_p_shuffl_modidx_areacurv_Sham>0)));
fig=figure(); hold on;
mod_idx_AreaCurveStim(find(mod_idx_AreaCurveStim<-1))=-0.5; %exceptional correction to change one data point
histogram(mod_idx_AreaCurveStim(find(All_p_shuffl_modidx_areacurv_Stim>0)),60,'DisplayStyle','stairs','Normalization','probability','EdgeColor','r')
histogram(mod_idx_AreaCurveSham(find(All_p_shuffl_modidx_areacurv_Sham>0)),60,'DisplayStyle','stairs','Normalization','probability','EdgeColor','b')
line(0*[1 1],ylim,'linestyle',':','color','k')
xlim([-1 1])
xlabel('Modulation Index (a.u.)')
ylabel('count #')
legend({'Stim','Sham'},'box','off')
title({['AreaCurve-Stim-vs-Sham'];['ranksum p=' num2str(p)]})
set(gca,'fontweight','bold','fontsize',11,'linewidth',1)
set(fig,'Renderer', 'painters')
saveas(gca,fullfile(cd,['CompareDistModIdx_AreaCurveStimVsSham.pdf']))
saveas(gca,fullfile(cd,['CompareDistModIdx_AreaCurveStimVsSham.svg']))
saveas(gca,fullfile(cd,['CompareDistModIdx_AreaCurveStimVsSham.fig']))
close all


