bachelor_thesis/DataAnalysisToolbox/reviewTimeEffects.m

%% Analyse potential time effects on the three main criteria
% This script offers multiple options to analyse difference of the first
% and last trial. The three main features interneuron spike amplitude, T
% cell respose to a -2 nA stimulus and the latency of the interneuron
% respose to a T cell spike can be analysed. 
% The first part offers curve sets for the T cell amplitudes and the spike 
% triggered average as well as trend curves for the interneuron spike 
% amplitude and the T cell amplitudes.
% The second part shows the 3D scatter plots of the three criteria. Each
% scatter plot is available three times: for the first trial, for the last
% trial and for the average over all trials.
% The third part offers a table that is printed to the command window
% giving the differences between the first and the last trial stated as
% numeric differences of the three criteria.
% All single appraoches can be used to estimate changes over the trials and
% to spot abnormalties, which then can be reviewed further. 
% -------------------------------------------------------------------------
% Program by:   Bjarne Schultze     [last modified: 10.03.2022]
% -------------------------------------------------------------------------

%% Analyse time effects with curve sets
% check the curvesets for the amplitudes (brighter colors indicate later
% trials)
ampAtStim_curveset(1:62, 'StimCell', 'INT')

% check the curvesets for the spike triggered average of the postsynaptic
% response
spikeTrigAvg_postsyn_curveset(41:62)

%% Analyse time effects with trend curves
% check the variability of the spike amplitude (and/or its measurement)
% for all trials together: measured spike height plotted against its index
intSpikeHeight_trend(1:62)

% create a trend curve for the measured amplitudes for the T cell during  
% the -2 and the 1.5 nA stimulus
% obtain the necessary data
amps = ampAtStim(1:62, 'Stimuli', [-2,1.5], 'StimCell', 'INT');

% plot the data in one figure per trial
for i = 1:62
    figure();
    stimINT = setdiff((1:size(amps{i}.AmplitudeT,2)),amps{i}.stimulatedT);
    xData = 1:length(stimINT);
    yData2nA = amps{i}.AmplitudeT(1,stimINT);
    yData15nA = amps{i}.AmplitudeT(2,stimINT);
    % create plot and labels
    plot(xData, yData2nA, 'o', 'Color', [0,0.34,0.64], ...
        'MarkerFaceColor',[0,0.34,0.64])
    hold on; plot(xData, yData15nA, 'o', 'Color', [1,0.7,0], ...
        'MarkerFaceColor', [1,0.7,0]); hold off
    xlabel("trial"); ylabel("T cell amplitude [mV]");
    title(sprintf("T cell response over time - File: %g",i))
    legend("response to -2 nA stimulus", "response to 1.5 nA stimulus")
end


%% Compare the three main features for the first and the last trial
% calculate the amplitudes for the T cell during int stimulation (-2 and
% 1.5 nA stimuli)
if ~exist('amps','var')
    amps = ampAtStim(1:62, 'Stimuli', [-2,1.5], 'StimCell', 'INT');
end
% pre-allocate vectors for data gathering
spikeHeights_firstTrial = zeros(62,1);
spikeHeights_lastTrial = zeros(62,1);
latencies_firstTrial = zeros(62,1);
latencies_lastTrial = zeros(62,1);
amp_m2nAStim_firstTrial = zeros(62,1);
amp_1_5nAStim_firstTrial = zeros(62,1);
amp_m2nAStim_lastTrial = zeros(62,1);
amp_1_5nAStim_lastTrial = zeros(62,1);

for i = 1:62
    % obtain the numbers of the first and the last trial
    stimT = amps{i}.stimulatedT;
    trials = 1:size(amps{i}.AmplitudeT,2);
    stimINT = setdiff(trials,stimT);
    firstTrialINT = min(stimINT);
    lastTrialINT = max(stimINT);
    firstTrialT = min(stimT);
    lastTrialT = max(stimT);
    % get the latencies for the int response in first and last trial
    latencies_firstTrial(i,1) = spikeTrigAvg_postsyn(i,firstTrialT);
    latencies_lastTrial(i,1) = spikeTrigAvg_postsyn(i, lastTrialT);
    % get the spike amplitudes for first and last trial
    spikeHeights_firstTrial(i,1) = intSpikeHeight(i,firstTrialINT);
    spikeHeights_lastTrial(i,1) = intSpikeHeight(i,lastTrialINT);
    % gather the T cell amplitudes for -2 and 1.5 nA stimuli
    amp_m2nAStim_firstTrial(i,1) = amps{i,1}.AmplitudeT(1,firstTrialINT);
    amp_1_5nAStim_firstTrial(i,1) = amps{i,1}.AmplitudeT(2,firstTrialINT);
    amp_m2nAStim_lastTrial(i,1) = amps{i,1}.AmplitudeT(1,lastTrialINT);
    amp_1_5nAStim_lastTrial(i,1) = amps{i,1}.AmplitudeT(2,lastTrialINT);
end

% correct the values for the non-spiking cells
spikeHeights_firstTrial([5,10,26,36,38,40]) = -10; 
spikeHeights_lastTrial([5,10,26,36,38,40]) = -10; 
latencies_lastTrial(isnan(latencies_lastTrial)) = 60;
latencies_firstTrial(isnan(latencies_firstTrial)) = 60;

% load additional files
load("../AdditionalFiles/decisionData.mat")
load("../AdditionalFiles/deanonymization.mat")

% get the numbers of the cells in each package
pack_1_3 = find(eq(deanonymization.package(:,1),13));
pack_2 = find(eq(deanonymization.package(:,1),2));
pack_4_6 = find(eq(deanonymization.package(:,1),46));
pack_5 = find(eq(deanonymization.package(:,1),5));

% take the assigned cell types into account
int218 = find(eq(deanonymization.cellType(:,1),218));
int162 = find(eq(deanonymization.cellType(:,1),162));
int264 = find(eq(deanonymization.cellType(:,1),264));
int212 = find(eq(deanonymization.cellType(:,1),212));
int201203 = find(eq(deanonymization.cellType(:,1),201203));
int159 = find(eq(deanonymization.cellType(:,1),159));
int157 = find(eq(deanonymization.cellType(:,1),157));
int6061 = find(eq(deanonymization.cellType(:,1),6061));

%% Analyse time effects by comparing the 3D scatter plots
% ---- scatter plots for all packages together ----------------------------
fig10 = figure();
tiledlayout(1,3)
fontSze = 11;

nexttile
scatter3(spike_heights(pack_1_3),amp_minus2nA_stim(pack_1_3), ...
    latencies(pack_1_3), 'MarkerEdgeColor',[0,0,0], 'MarkerFaceColor', ...
    [0.2,0.6,1])
hold on;
scatter3(spike_heights(pack_2),amp_minus2nA_stim(pack_2),latencies(pack_2), ...
    'MarkerEdgeColor',[0,0,0], 'MarkerFaceColor',[1,0.6,0.2])
scatter3(spike_heights(pack_4_6),amp_minus2nA_stim(pack_4_6), ...
    latencies(pack_4_6), 'MarkerEdgeColor',[0,0,0], 'MarkerFaceColor', ...
    [0.6,0.8,0.2])
scatter3(spike_heights(pack_5),amp_minus2nA_stim(pack_5),latencies(pack_5), ...
    'MarkerEdgeColor',[0,0,0], 'MarkerFaceColor',[0.6,0.2,0.6])
hold off;
% labelling
xlabel("spike amplitude [mV]");
ylabel("T cell amplitude [mV]");
zlabel("latency [ms]");
title("all trials")
legend("Packages 1 & 3", "Package 2", "Packages 4 & 6", ...
    "Packages 5", 'FontSize', fontSze);

nexttile
scatter3(spikeHeights_firstTrial(pack_1_3),...
    amp_m2nAStim_firstTrial(pack_1_3), latencies_firstTrial(pack_1_3), ...
    'MarkerEdgeColor',[0,0,0], 'MarkerFaceColor', ...
    [0.2,0.6,1])
hold on;
scatter3(spikeHeights_firstTrial(pack_2), ...
    amp_m2nAStim_firstTrial(pack_2),latencies_firstTrial(pack_2), ...
    'MarkerEdgeColor',[0,0,0], 'MarkerFaceColor',[1,0.6,0.2])
scatter3(spikeHeights_firstTrial(pack_4_6),amp_minus2nA_stim(pack_4_6), ...
    latencies_firstTrial(pack_4_6), 'MarkerEdgeColor',[0,0,0], ...
    'MarkerFaceColor', ...
    [0.6,0.8,0.2])
scatter3(spikeHeights_firstTrial(pack_5), ...
    amp_m2nAStim_firstTrial(pack_5),latencies_firstTrial(pack_5), ...
    'MarkerEdgeColor',[0,0,0], 'MarkerFaceColor',[0.6,0.2,0.6])
hold off;
% labelling
xlabel("spike amplitude [mV]");
ylabel("T cell amplitude [mV]");
zlabel("latency [ms]");
title("first trial")

nexttile
scatter3(spikeHeights_lastTrial(pack_1_3), ...
    amp_m2nAStim_lastTrial(pack_1_3), latencies_lastTrial(pack_1_3), ...
    'MarkerEdgeColor',[0,0,0], 'MarkerFaceColor', ...
    [0.2,0.6,1])
hold on;
scatter3(spikeHeights_lastTrial(pack_2), ...
    amp_m2nAStim_lastTrial(pack_2),latencies_lastTrial(pack_2), ...
    'MarkerEdgeColor',[0,0,0], 'MarkerFaceColor',[1,0.6,0.2])
scatter3(spikeHeights_lastTrial(pack_4_6), ...
    amp_m2nAStim_lastTrial(pack_4_6), latencies_lastTrial(pack_4_6), ...
    'MarkerEdgeColor',[0,0,0], 'MarkerFaceColor', ...
    [0.6,0.8,0.2])
scatter3(spikeHeights_lastTrial(pack_5), ...
    amp_m2nAStim_lastTrial(pack_5),latencies_lastTrial(pack_5), ...
    'MarkerEdgeColor',[0,0,0], 'MarkerFaceColor',[0.6,0.2,0.6])
hold off;
% labelling
xlabel("spike amplitude [mV]");
ylabel("T cell amplitude [mV]");
zlabel("latency [ms]");
title("last trial")

% ---- scatter plots for packages 4 to 5 ----------------------------------
% set up the figure
fig11 = figure;
tiledlayout(2,3)

% plot the data
nexttile                            % packages 4 and 6 (interneurons 60-61)
scatter3(spike_heights(int6061),amp_minus2nA_stim(int6061),...
    latencies(int6061), 'MarkerEdgeColor',[0,0,0], ...
    'MarkerFaceColor',[0.87,0.87,0.23]);
title("all trials")
% labelling
xlabel("spike amplitude [mV]");
ylabel("T cell amplitude [mV]");
zlabel("latency [ms]");
legend("int 60-61")

nexttile                            % packages 4 and 6 (interneurons 60-61)
scatter3(spikeHeights_firstTrial(int6061),amp_m2nAStim_firstTrial(int6061),...
    latencies_firstTrial(int6061), 'MarkerEdgeColor',[0,0,0], ...
    'MarkerFaceColor',[0.87,0.87,0.23]);
title("first trial")
% labelling
xlabel("spike amplitude [mV]");
ylabel("T cell amplitude [mV]");
zlabel("latency [ms]");
legend("int 60-61")

nexttile
scatter3(spikeHeights_lastTrial(int6061),amp_m2nAStim_lastTrial(int6061),...
    latencies_lastTrial(int6061), 'MarkerEdgeColor',[0,0,0], ...
    'MarkerFaceColor',[0.87,0.87,0.23]);
title("last trial")
% labelling
xlabel("spike amplitude [mV]");
ylabel("T cell amplitude [mV]");
zlabel("latency [ms]");
legend("int 60-61")

nexttile                            % package 5 (interneuron 264)
scatter3(spike_heights(int264),amp_minus2nA_stim(int264),...
    latencies(int264), 'MarkerEdgeColor',[0,0,0], ...
    'MarkerFaceColor',[0,0.4,0.7]);
title("all trials")
% labelling
xlabel("spike amplitude [mV]");
ylabel("T cell amplitude [mV]");
zlabel("latency [ms]");
legend("int 264")

nexttile                             % package 5 (interneuron 264)
scatter3(spikeHeights_firstTrial(int264),amp_m2nAStim_firstTrial(int264),...
    latencies_firstTrial(int264),'MarkerEdgeColor',[0,0,0], ...
    'MarkerFaceColor',[0,0.4,0.7]);
title("first trial")
% labelling
xlabel("spike amplitude [mV]");
ylabel("T cell amplitude [mV]");
zlabel("latency [ms]");
legend("int 264")

nexttile                             % package 5 (interneuron 264)
scatter3(spikeHeights_lastTrial(int264),amp_m2nAStim_lastTrial(int264),...
    latencies_lastTrial(int264),'MarkerEdgeColor',[0,0,0], ...
    'MarkerFaceColor',[0,0.4,0.7]);
title("last trial")
% labelling
xlabel("spike amplitude [mV]");
ylabel("T cell amplitude [mV]");
zlabel("latency [ms]");
legend("int 264")

% adjust axes
ax = fig11.Children.Children;
ax(2).XLim = [-10,60]; ax(2).YLim = [-10,0]; ax(2).ZLim = [0,60];
ax(4).XLim = [-10,60]; ax(4).YLim = [-10,0]; ax(4).ZLim = [0,60];
ax(6).XLim = [-10,60]; ax(6).YLim = [-10,0]; ax(6).ZLim = [0,60];
ax(8).XLim = [-10,60]; ax(8).YLim = [-10,0]; ax(8).ZLim = [0,60];
ax(10).XLim = [-10,60]; ax(12).YLim = [-10,0]; ax(12).ZLim = [0,60];
ax(12).XLim = [-10,60]; ax(12).YLim = [-10,0]; ax(12).ZLim = [0,60];

% ---- scatter plot for packge 2 ------------------------------------------
% plot the data
fig12 = figure;
tiledlayout(1,3)

nexttile
scatter3(spike_heights(int218), amp_minus2nA_stim(int218), ...
    latencies(int218),'MarkerEdgeColor',[0,0,0], ...
    'MarkerFaceColor',[0.4,0.8,0.8])
hold on;
scatter3(spike_heights(int212), amp_minus2nA_stim(int212), ...
    latencies(int212),'MarkerEdgeColor',[0,0,0], ...
    'MarkerFaceColor',[0.8,0.2,0.1])
scatter3(spike_heights(int201203),amp_minus2nA_stim(int201203),...
    latencies(int201203), 'MarkerEdgeColor',[0,0,0],...
    'MarkerFaceColor', [0.4,0.7,0.1])
hold off;
% labelling
xlabel("spike amplitude [mV]");
ylabel("T cell amplitude [mV]");
zlabel("latency [ms]");
legend("int 218","int 212", "int 201-203")
title("all trials")

nexttile
scatter3(spikeHeights_firstTrial(int218), amp_m2nAStim_firstTrial(int218), ...
    latencies_firstTrial(int218),'MarkerEdgeColor',[0,0,0], ...
    'MarkerFaceColor',[0.4,0.8,0.8])
hold on;
scatter3(spikeHeights_firstTrial(int212), amp_m2nAStim_firstTrial(int212), ...
    latencies_firstTrial(int212),'MarkerEdgeColor',[0,0,0], ...
    'MarkerFaceColor',[0.8,0.2,0.1])
scatter3(spikeHeights_firstTrial(int201203), ...
    amp_m2nAStim_firstTrial(int201203), latencies_firstTrial(int201203), ...
    'MarkerEdgeColor',[0,0,0],'MarkerFaceColor', [0.4,0.7,0.1])
hold off;
% labelling
xlabel("spike amplitude [mV]");
ylabel("T cell amplitude [mV]");
zlabel("latency [ms]");
title("first trial")

nexttile
scatter3(spikeHeights_lastTrial(int218), amp_m2nAStim_lastTrial(int218), ...
    latencies_lastTrial(int218),'MarkerEdgeColor',[0,0,0], ...
    'MarkerFaceColor',[0.4,0.8,0.8])
hold on;
scatter3(spikeHeights_lastTrial(int212), amp_m2nAStim_lastTrial(int212), ...
    latencies_lastTrial(int212),'MarkerEdgeColor',[0,0,0], ...
    'MarkerFaceColor',[0.8,0.2,0.1])
scatter3(spikeHeights_lastTrial(int201203), ...
    amp_m2nAStim_lastTrial(int201203), latencies_lastTrial(int201203), ...
    'MarkerEdgeColor',[0,0,0],'MarkerFaceColor', [0.4,0.7,0.1])
hold off;
% labelling
xlabel("spike amplitude [mV]");
ylabel("T cell amplitude [mV]");
zlabel("latency [ms]");
title("last trial")

% adjust axes
ax = fig12.Children.Children;

ax(1).XLim = [-10,60]; ax(1).YLim = [-10,0]; ax(1).ZLim = [0,60];
ax(2).XLim = [-10,60]; ax(2).YLim = [-10,0]; ax(2).ZLim = [0,60];
ax(4).XLim = [-10,60]; ax(4).YLim = [-10,0]; ax(4).ZLim = [0,60];


% ---- scatter plots for packages 1 and 3 ---------------------------------
% plot the data
fig13 = figure;
tiledlayout(1,3)

nexttile
scatter3(spike_heights(int159), amp_minus2nA_stim(int159), ...
    latencies(int159),'MarkerEdgeColor',[0,0,0], ...
    'MarkerFaceColor',[0.1,0.8,0.2])
hold on;
scatter3(spike_heights(int157), amp_minus2nA_stim(int157), ...
    latencies(int157),'MarkerEdgeColor',[0,0,0], ...
    'MarkerFaceColor',[1,0.75,0])
scatter3(spike_heights(int162), amp_minus2nA_stim(int162), ...
    latencies(int162),'MarkerEdgeColor',[0,0,0], ...
    'MarkerFaceColor',[0.7,0.7,0.4])
hold off;
% labelling
xlabel("spike amplitude [mV]");
ylabel("T cell amplitude [mV]");
zlabel("latency [ms]");
legend("int 159","int 157", "int 162")

nexttile
scatter3(spikeHeights_firstTrial(int159), amp_m2nAStim_firstTrial(int159), ...
    latencies_firstTrial(int159),'MarkerEdgeColor',[0,0,0], ...
    'MarkerFaceColor',[0.1,0.8,0.2])
hold on;
scatter3(spikeHeights_firstTrial(int157), amp_m2nAStim_firstTrial(int157), ...
    latencies_firstTrial(int157),'MarkerEdgeColor',[0,0,0], ...
    'MarkerFaceColor',[1,0.75,0])
scatter3(spikeHeights_firstTrial(int162), amp_m2nAStim_firstTrial(int162), ...
    latencies_firstTrial(int162),'MarkerEdgeColor',[0,0,0], ...
    'MarkerFaceColor',[0.7,0.7,0.4])
hold off;
% labelling
xlabel("spike amplitude [mV]");
ylabel("T cell amplitude [mV]");
zlabel("latency [ms]");
title("first trial")

nexttile
scatter3(spikeHeights_lastTrial(int159), amp_m2nAStim_lastTrial(int159), ...
    latencies_lastTrial(int159),'MarkerEdgeColor',[0,0,0], ...
    'MarkerFaceColor',[0.1,0.8,0.2])
hold on;
scatter3(spikeHeights_lastTrial(int157), amp_m2nAStim_lastTrial(int157), ...
    latencies_lastTrial(int157),'MarkerEdgeColor',[0,0,0], ...
    'MarkerFaceColor',[1,0.75,0])
scatter3(spikeHeights_lastTrial(int162), amp_m2nAStim_lastTrial(int162), ...
    latencies_lastTrial(int162),'MarkerEdgeColor',[0,0,0], ...
    'MarkerFaceColor',[0.7,0.7,0.4])
hold off;
% labelling
xlabel("spike amplitude [mV]");
ylabel("T cell amplitude [mV]");
zlabel("latency [ms]");
title("last trial")

% adjust the axes
ax = fig13.Children.Children;
ax(1).XLim = [-10,60]; ax(1).YLim = [-10,0]; ax(1).ZLim = [0,60];
ax(2).XLim = [-10,60]; ax(2).YLim = [-10,0]; ax(2).ZLim = [0,60];
ax(4).XLim = [-10,60]; ax(4).YLim = [-10,0]; ax(4).ZLim = [0,60];

%% Analyse time effect by the variation of the values, comparison table
% calculate the differences for first and last trial
diff_spikeAmps = abs(spikeHeights_firstTrial-spikeHeights_lastTrial);
diff_amp2nAStim = abs(amp_m2nAStim_firstTrial-amp_m2nAStim_lastTrial);
diff_latencies = abs(latencies_firstTrial-latencies_lastTrial);

% create a header for the output
feature_list(1,1) = sprintf("\nComparison of the main features " + ...
    "regarding their difference from first to last trial\n" + ...
    "---------------------------------------------\n" + ...
    "Cell \t\x0394 spike \t\x0394 amp -2 nA \t\x0394 latency\n" + ...
    "number \tamplitude\tstimulus \n" + ...
    "\t\t[mV] \t\t[mV] \t\t\t[ms] \n" + ...
    "---------------------------------------------\n");
% write the difference values into seperate lines with a nice layout
for i = 1:62
    feature_list(i+1,1) = sprintf("%02.f | \t%05.2f \t\t%04.2f " + ...
        "\t\t\t%04.1f\n", ...
        i, diff_spikeAmps(i,1), diff_amp2nAStim(i,1), ...
        diff_latencies(i,1));
end

% add an endline for the output table
feature_list(64,1) = ...
    sprintf("---------------------------------------------\n");
% print the table to the command window
fprintf("%s", feature_list)