load('PooledCellsData.mat')
%% Connected dots - Hit trials
figure

subplot(121)
ConnectedDots2 (Control_hit_freq(:,1), Control_hit_freq(:,2)) % Unclear_hit_freq(:,1)])
%set(gca,'XTickLabel',{'Control (n=41)','hM4Di (n=40)'})
ylabel('Spike rate (Hz)')
title('Control')
[p h] = signrank (Control_hit_freq(:,1), Control_hit_freq(:,2)) % Unclear_hit_freq(:,1)])
hold on
text(1.2, 20, ['p = ', num2str(p)]) 

subplot(122)
ConnectedDots2 (newDREADD_hit_freq(:,1), newDREADD_hit_freq(:,2)) % Unclear_hit_freq(:,1)])
%set(gca,'XTickLabel',{'Control (n=42)','hM4Di (n=26)'})
ylabel('Spike rate (Hz)')
title('hM4Di')
[p h] = signrank (newDREADD_hit_freq(:,1), newDREADD_hit_freq(:,2)) % Unclear_hit_freq(:,1)])
hold on
text(1.2, 20, ['p = ', num2str(p)])  
ylim([0 25])

figure

subplot(121)
ConnectedDots2 (Control_hit_burst(:,1), Control_hit_burst(:,2)) % Unclear_hit_burst(:,1)])
%set(gca,'XTickLabel',{'Control (n=42)','hM4Di (n=26)'})
ylabel('Burst rate (Hz)')
title('Control')
[p h] = signrank (Control_hit_burst(:,1), Control_hit_burst(:,2)) % Unclear_hit_burst(:,1)])
hold on
text(1.2, 2, ['p = ', num2str(p)]) 

subplot(122)
ConnectedDots2 (newDREADD_hit_burst(:,1), newDREADD_hit_burst(:,2)) % Unclear_hit_burst(:,1)])
%set(gca,'XTickLabel',{'Control (n=42)','hM4Di (n=26)'})
ylabel('Burst rate (Hz)')
title('hM4Di')
[p h] = signrank (newDREADD_hit_burst(:,1), newDREADD_hit_burst(:,2)) % Unclear_hit_burst(:,1)])
hold on
text(1.2, 2, ['p = ', num2str(p)])  
ylim([0 2.5])

%% Connected Dots - miss trials
figure

subplot(121)
ConnectedDots2 (Control_miss_freq(:,1), Control_miss_freq(:,2)) % Unclear_miss_freq(:,1)])
%set(gca,'XTickLabel',{'Control (n=41)','hM4Di (n=40)'})
ylabel('Spike rate (Hz)')
title('Control')
ylim ([0 25])
[p h] = signrank (Control_miss_freq(:,1), Control_miss_freq(:,2)) % Unclear_miss_freq(:,1)])
hold on
text(1.2, 20, ['p = ', num2str(p)]) 

subplot(122)
ConnectedDots2 (newDREADD_miss_freq(:,1), newDREADD_miss_freq(:,2)) % Unclear_miss_freq(:,1)])
%set(gca,'XTickLabel',{'Control (n=42)','hM4Di (n=26)'})
ylabel('Spike rate (Hz)')
title('hM4Di')
[p h] = signrank (newDREADD_miss_freq(:,1), newDREADD_miss_freq(:,2)) % Unclear_miss_freq(:,1)])
hold on
text(1.2, 20, ['p = ', num2str(p)])  
ylim([0 25])

figure

subplot(121)
ConnectedDots2 (Control_miss_burst(:,1), Control_miss_burst(:,2)) % Unclear_miss_burst(:,1)])
%set(gca,'XTickLabel',{'Control (n=42)','hM4Di (n=26)'})
ylabel('Burst rate (Hz)')
title('Control')
ylim([0 2.5])
[p h] = signrank (Control_miss_burst(:,1), Control_miss_burst(:,2)) % Unclear_miss_burst(:,1)])
hold on
text(1.2, 2, ['p = ', num2str(p)]) 

subplot(122)
ConnectedDots2 (newDREADD_miss_burst(:,1), newDREADD_miss_burst(:,2)) % Unclear_miss_burst(:,1)])
%set(gca,'XTickLabel',{'Control (n=42)','hM4Di (n=26)'})
ylabel('Burst rate (Hz)')
title('hM4Di')
[p h] = signrank (newDREADD_miss_burst(:,1), newDREADD_miss_burst(:,2)) % Unclear_miss_burst(:,1)])
hold on
text(1.2, 2, ['p = ', num2str(p)])  
ylim([0 2.5])


%% Difference

Control_hit_freq(:,3) = Control_hit_freq(:,2) - Control_hit_freq(:,1);
newDREADD_hit_freq(:,3) = newDREADD_hit_freq(:,2) - newDREADD_hit_freq(:,1);
Unclear_hit_freq(:,3) =Unclear_hit_freq(:,2) - Unclear_hit_freq(:,1);

figure
subplot(121)
PlotAllDataPoints (Control_hit_freq(:,3), newDREADD_hit_freq(:,3) )
set(gca,'XTickLabel',{'Control (n=42)','hM4Di (n=26)'})
ylabel('Post - Pre')
title('Spike rate change')

[p h] = ranksum (Control_hit_freq(:,3), newDREADD_hit_freq(:,3))

hold on
text(1.2, 10, ['p = ', num2str(p)]) 

subplot(122)
PlotAllDataPoints (abs(Control_hit_freq(:,3)), abs(newDREADD_hit_freq(:,3)) )
set(gca,'XTickLabel',{'Control (n=42)','hM4Di (n=26)'})
ylabel('|Post - Pre|')
title('Spike rate Absolute change')

[p h] = ranksum (abs(Control_hit_freq(:,3)), abs(newDREADD_hit_freq(:,3)))

hold on
text(1.2, 10, ['p = ', num2str(p)])

% Bursts
Control_hit_burst(:,3) = Control_hit_burst(:,2) - Control_hit_burst(:,1);
newDREADD_hit_burst(:,3) = newDREADD_hit_burst(:,2) - newDREADD_hit_burst(:,1);
Unclear_hit_burst(:,3) =Unclear_hit_burst(:,2) - Unclear_hit_burst(:,1);


figure
subplot(121)
PlotAllDataPoints (Control_hit_burst(:,3), newDREADD_hit_burst(:,3) )
set(gca,'XTickLabel',{'Control (n=42)','hM4Di (n=26)'})
ylabel('Post - Pre')
title('Burst rate change')

[p h] = ranksum (Control_hit_burst(:,3), newDREADD_hit_burst(:,3))

hold on
text(1.2, 1, ['p = ', num2str(p)]) 


subplot(122)
PlotAllDataPoints (abs(Control_hit_burst(:,3)), abs(newDREADD_hit_burst(:,3)) )
set(gca,'XTickLabel',{'Control (n=42)','hM4Di (n=26)'})
ylabel('|Post - Pre|')
title('Spike rate Absolute change')

[p h] = ranksum (abs(Control_hit_burst(:,3)), abs(newDREADD_hit_burst(:,3)))

hold on
text(1.2, 1, ['p = ', num2str(p)]) 
%% Relative change
Control_baseline = mean(Control_hit_burst(:,1));
newDREADD_baseline = mean(newDREADD_hit_burst(:,1));

Control_hit_burst(:,4) = Control_hit_burst(:,3)./Control_baseline;
newDREADD_hit_burst(:,4) = newDREADD_hit_burst(:,3)./newDREADD_baseline;

[p h] = ranksum (Control_hit_burst(:,4), newDREADD_hit_burst(:,4))

figure
PlotAllDataPoints (Control_hit_burst(:,4).*100, newDREADD_hit_burst(:,4).*100 )
set(gca,'XTickLabel',{'Control (n=42)','hM4Di (n=26)'})
ylabel('Relative change (%)')
title('Burst')


%% Post stim
figure

subplot(121)
PlotAllDataPoints (Control_hit_freq(:,2), newDREADD_hit_freq(:,2) )
set(gca,'XTickLabel',{'Control (n=41)','hM4Di (n=41)'})
ylabel('Spike rate (Hz)')
title('Hit - Post stimulation')

[p h] = ranksum (Control_hit_freq(:,2), newDREADD_hit_freq(:,2))

hold on
text(1.2, 10, ['p = ', num2str(p)]) 

subplot(122)
PlotAllDataPoints (Control_miss_freq(:,2), newDREADD_miss_freq(:,2) )
set(gca,'XTickLabel',{'Control (n=18)','hM4Di (n=17)'})
ylabel('Spike rate (Hz)')
title('Miss - Post stimulation')

[p h] = ranksum (Control_miss_freq(:,2), newDREADD_miss_freq(:,2))

hold on
text(1.2, 8, ['p = ', num2str(p)]) 
ylim([0 25])

%% Bar wit error
N1 = size(Control_hit_freq,1)
Control_Avg_hitFreq = nanmean(Control_hit_freq);
Control_SEM_hitFreq = nanstd(Control_hit_freq,1)./sqrt(N1);
Control_Avg_hitBurst = nanmean(Control_hit_burst)  ;
Control_SEM_hitBurst = nanstd(Control_hit_burst,1)./sqrt(N1);

N2 = size(newDREADD_hit_freq,1)
newDREADD_Avg_hitFreq = nanmean(newDREADD_hit_freq);
newDREADD_SEM_hitFreq = nanstd(newDREADD_hit_freq,1)./sqrt(N2);
newDREADD_Avg_hitBurst = nanmean(newDREADD_hit_burst)  ;
newDREADD_SEM_hitBurst = nanstd(newDREADD_hit_burst,1)./sqrt(N2);


% %% Plotting
figure
subplot(2,1,1)
varargout = barwitherr([Control_SEM_hitFreq; newDREADD_SEM_hitFreq],[Control_Avg_hitFreq; newDREADD_Avg_hitFreq])  
set(gca,'XTickLabel',{'Control (n=41)','hM4Di (n=40)'})
ylabel('Firing rate (Hz)')
ylim([0 5])

title('hit trials')
legend({'preStim', 'postStim'})

subplot(2,1,2)
varargout = barwitherr([Control_SEM_hitBurst; newDREADD_SEM_hitBurst],[Control_Avg_hitBurst;newDREADD_Avg_hitBurst])  
set(gca,'XTickLabel',{'Control (n=41)','hM4Di (n=40)'})
ylabel('Burst rate (Hz)')


%% Violin plot
figure
violin(Control_hit_burst);
set(gca, 'XTick',[1,2],'XTickLabel',{'PreStim', 'PostStim'})
title('Control')

figure
violin(newDREADD_hit_burst);
set(gca, 'XTick',[1,2],'XTickLabel',{'PreStim', 'PostStim'})
title('hM4Di')