VentralPallidumAlcohol2019/AnalysisAndFigures/a_BehaviorInstumentalTrainingAnalysis.m

%% Import data from excel
clear;

[~,~,PAInstTrainingData] = xlsread('PAInstTrainingDataALL.xlsx');

%% Assign data to corresponding variables within a structure
VarNames = PAInstTrainingData(1,:);
Data = PAInstTrainingData(2:end,:);

PAInst = struct();

for i=1:42
    PAInst.(VarNames{i}) = Data(1:end,(i));
end

PAInst.CSPlusOnset = Data(1:end,43:103);
PAInst.CSMinOnset = Data(1:end,104:163);
PAInst.PETimestamps = Data(1:end,225:942);

%% Convert rat numbers from strings into numbers and append to a separate field
for i = 1:length(PAInst.Rat)
   PAInst(:).RatNum(i,1) = str2num(PAInst.Rat{i}(3:5));
end

%% Sort matrix by rat
[sortedRat,index] = sort([PAInst(:).Rat]);
fields = fieldnames(PAInst);
for i = 1:length(fields)
    PAInst.(fields{i}) = PAInst.(fields{i})(index,:);
end

%% Number the start dates for each rat and display any errors
Rat = PAInst.Rat(1);
k = 1;
for i=1:length(PAInst.StartDate)
   if  strcmp(PAInst.Rat(i), Rat)
       PAInst.Day(i,:) = k;
       k=k+1;
       if i < length(PAInst.StartDate)
           if PAInst.StartDate{i} == PAInst.StartDate{i+1}
               fprintf('repeated date %d with rat %s\n', PAInst.StartDate{i}, PAInst.Rat{i});
           end
       end
   else
       k = 1;
       Rat = PAInst.Rat(i);
       PAInst.Day(i,:) = k;
       k=k+1;
   end
end

%% Determine stages of training
% look at durations
% if duration = 60 ==> stage 1
% if duration = 30 ==> stage 2
% if duration = 20 ==> stage 3
% if duration = 10 and numtrials <= 60 ==> stage 4
% if duration = 10 and numtrials > 60 ==> stage 5
% if duration = anything else ==> stage NaN, print rat and date in which this occurs

for i=1:length(PAInst.Duration)
    if PAInst.Duration{i} >= 60
        PAInst.Stage{i,:} = 1;
    elseif (PAInst.Duration{i} == 30 || PAInst.Duration{i} == 40)
        PAInst.Stage{i,:} = 2;
    elseif PAInst.Duration{i} == 20
        PAInst.Stage{i,:} = 3;
    elseif PAInst.Duration{i} == 10 & PAInst.CSMinOnset{i,1}==0
        PAInst.Stage{i,:} = 4;
    elseif PAInst.Duration{i} == 10 & PAInst.CSMinOnset{i,1}>0
        PAInst.Stage{i,:} = 5;
    else
        PAInst.Stage{i,:} = NaN;
        fprintf('stage is NaN for rat %s on date %d\n', PAInst.Rat{i}, PAInst.StartDate{i});   
    end    
end

%% Setting NS values to 0

for i=1:length(PAInst.CSMinRatio)
    if PAInst.Stage{i} < 5
        PAInst.CSMinRatio{i} = NaN;
    end
end

%% Calculate trial-by-trial PE latencies for CSPlus and CSMin

% Convert necessary cells to matrices, setting all strings to NaN
for y=1:length(PAInst.CSPlusOnset)
    for x=1:size(PAInst.CSPlusOnset,2)
        if ischar(PAInst.CSPlusOnset{y,x})
           PAInst.CSPlusOnset{y,x} = NaN;
        end
    end
end

for y=1:length(PAInst.CSMinOnset)
    for x=1:size(PAInst.CSMinOnset,2)
        if ischar(PAInst.CSMinOnset{y,x})
           PAInst.CSMinOnset{y,x} = NaN;
        end
    end
end

for y=1:size(PAInst.PETimestamps,1)
    for x=1:size(PAInst.PETimestamps,2)
        if ischar(PAInst.PETimestamps{y,x})
           PAInst.PETimestamps{y,x} = NaN;
        end
    end
end

PAInst.CSPlusOnset = cell2mat(PAInst.CSPlusOnset);
PAInst.CSMinOnset = cell2mat(PAInst.CSMinOnset);
PAInst.PETimestamps = cell2mat(PAInst.PETimestamps);

% Calculate CSPlus absolute and cue period latencies
%   Conditions:
%       if next csplus > 0 but less than pe time, lat = nan
%       if next csmin > 0 but less than pe time, lat =nan

for y=1:length(PAInst.CSPlusOnset)
    for x=1:size(PAInst.CSPlusOnset,2)
       nextPlus = find(PAInst.CSPlusOnset(y,:) > PAInst.CSPlusOnset(y,x), 1);
       nextMin = find(PAInst.CSMinOnset(y,:) > PAInst.CSPlusOnset(y,x), 1);
       firstPE = find(PAInst.PETimestamps(y,:) > PAInst.CSPlusOnset(y,x), 1);
       PAInst.CSPlusAbsLat{y,x} = PAInst.PETimestamps(y,firstPE) - PAInst.CSPlusOnset(y,x);
       if PAInst.CSPlusAbsLat{y,x} < PAInst.Duration{y,1}
           PAInst.CSPlusCPLat{y,x} = PAInst.CSPlusAbsLat{y,x};
       else
           PAInst.CSPlusCPLat{y,x} = NaN;
       end
       if PAInst.CSMinOnset(y,nextMin) > 0 & PAInst.CSMinOnset(y,nextMin) < PAInst.PETimestamps(y,firstPE)
           PAInst.CSPlusAbsLat{y,x} = NaN;
           PAInst.CSPlusCPLat{y,x} = NaN;
       elseif PAInst.CSPlusOnset(y,nextPlus) > 0 & PAInst.CSPlusOnset(y,nextPlus) < PAInst.PETimestamps(y,firstPE)
           PAInst.CSPlusAbsLat{y,x} = NaN;
           PAInst.CSPlusCPLat{y,x} = NaN;
       elseif isnan(PAInst.CSPlusOnset(y,x))
           PAInst.CSPlusAbsLat{y,x} = NaN;
           PAInst.CSPlusCPLat{y,x} = NaN;
       elseif PAInst.CSPlusOnset(y,x) == 0
           PAInst.CSPlusAbsLat{y,x} = NaN;
           PAInst.CSPlusCPLat{y,x} = NaN;
       elseif isempty(firstPE)
           PAInst.CSPlusAbsLat{y,x} = NaN;
           PAInst.CSPlusCPLat{y,x} = NaN;
       end
    end
end

% Calculate CSPlus absolute and cue period latency means
PAInst.CSPlusAbsLat = cell2mat(PAInst.CSPlusAbsLat);
PAInst.CSPlusCPLat = cell2mat(PAInst.CSPlusCPLat);
PAInst.CSPlusAbsLatMean = nanmean(PAInst.CSPlusAbsLat,2);
PAInst.CSPlusCPLatMean = nanmean(PAInst.CSPlusCPLat,2);

% Calculate CSMin absolute and cue period latencies
%   Conditions:
%       if next csplus > 0 but less than pe time, lat = nan
%       if next csmin > 0 but less than pe time, lat =nan

for y=1:length(PAInst.CSMinOnset)
    for x=1:size(PAInst.CSMinOnset,2)
       nextPlus = find(PAInst.CSPlusOnset(y,:) > PAInst.CSMinOnset(y,x), 1);
       nextMin = find(PAInst.CSMinOnset(y,:) > PAInst.CSMinOnset(y,x), 1);
       firstPE = find(PAInst.PETimestamps(y,:) > PAInst.CSMinOnset(y,x), 1);
       PAInst.CSMinAbsLat{y,x} = PAInst.PETimestamps(y,firstPE) - PAInst.CSMinOnset(y,x);
       if PAInst.CSMinAbsLat{y,x} < PAInst.Duration{y,1}
           PAInst.CSMinCPLat{y,x} = PAInst.CSMinAbsLat{y,x};
       else
           PAInst.CSMinCPLat{y,x} = NaN;
       end
       if PAInst.CSMinOnset(y,nextMin) > 0 & PAInst.CSMinOnset(y,nextMin) < PAInst.PETimestamps(y,firstPE)
           PAInst.CSMinAbsLat{y,x} = NaN;
           PAInst.CSMinCPLat{y,x} = NaN;
       elseif PAInst.CSPlusOnset(y,nextPlus) > 0 & PAInst.CSPlusOnset(y,nextPlus) < PAInst.PETimestamps(y,firstPE)
           PAInst.CSMinAbsLat{y,x} = NaN;
           PAInst.CSMinCPLat{y,x} = NaN;
       elseif isnan(PAInst.CSMinOnset(y,x))
           PAInst.CSMinAbsLat{y,x} = NaN;
           PAInst.CSMinCPLat{y,x} = NaN;
       elseif PAInst.CSMinOnset(y,x) == 0
           PAInst.CSMinAbsLat{y,x} = NaN;
           PAInst.CSMinCPLat{y,x} = NaN;
       elseif isempty(firstPE)
           PAInst.CSMinAbsLat{y,x} = NaN;
           PAInst.CSMinCPLat{y,x} = NaN;
       end
    end
end

% Calculate CSMin absolute and cue period latency means
PAInst.CSMinAbsLat = cell2mat(PAInst.CSMinAbsLat);
PAInst.CSMinCPLat = cell2mat(PAInst.CSMinCPLat);
PAInst.CSMinAbsLatMean = nanmean(PAInst.CSMinAbsLat,2);
PAInst.CSMinCPLatMean = nanmean(PAInst.CSMinCPLat,2);


%% Determine when rats reach criteria
% FOR EACH RAT:
% criteria = 0 if cs+ratio < .7 or cs-ratio > .3
% if cs+ratio > .7 and cs-ratio < .3, then
%   for the first time this happens, criteria = 0
%   for every other time this happens and even if they don't meet criteria in later dates (don't need to check cs ratios), criteria = 1
% if new rat, start over

k = 1;
for i=1:length(PAInst.Stage)
    if PAInst.Stage{i} < 5
        criteria_met = false; 
        PAInst.Criteria{i,1} = 0;
    elseif isnan(PAInst.Stage{i})
        criteria_met = false;
        PAInst.Criteria{i,1} = NaN;
    elseif criteria_met == true % this condition is set from criteria being met
        PAInst.Criteria{i,1} = 1;
    elseif PAInst.Stage{i,1} == 5
        if PAInst.CSPlusRatio{i} >= 0.7 & PAInst.CSMinRatio{i} <= 0.3 % correct criteria
            PAInst.Criteria{i,1} = 0;
            PAInst.CriteriaDay(k,1) = PAInst.Day(i,1); 
            PAInst.CriteriaDayRat(k,1) = PAInst.RatNum(i,1);
            PAInst.CriteriaDayCSPlusRatio(k,1) = PAInst.CSPlusRatio(i,1);
            PAInst.CriteriaDayCSMinRatio(k,1) = PAInst.CSMinRatio(i,1);
            PAInst.CriteriaDayCSPlusLat(k,1) = PAInst.CSPlusCPLatMean(i,1);
            PAInst.CriteriaDayCSMinLat(k,1) = PAInst.CSMinCPLatMean(i,1);
            k = k+1;
            criteria_met = true; % false switched to true 
        else
            PAInst.Criteria{i,1} = 0;
            criteria_met = false;
        end
    end
end

PAInst.Criteria = cell2mat(PAInst.Criteria);
PAInst.CSPlusRatio = cell2mat(PAInst.CSPlusRatio);
PAInst.CSMinRatio = cell2mat(PAInst.CSMinRatio);

%% Calculate averages of latency means on days where criteria = 0 during cue period (group data)

% To find averages of CS+ and CS- Cue Period Latency Means by day

RatSelection=ismember(PAInst.RatNum,[105 107 139 151]);

for i=1:max(PAInst.Day(PAInst.Criteria == 0))
    PAInst.AverageCSPlusCPLatMean(i) = nanmean(PAInst.CSPlusCPLatMean(PAInst.Criteria == 0 & PAInst.Day == i & RatSelection));
    PAInst.Error_AverageCSPlusCPLatMean(i) = nanste(PAInst.CSPlusCPLatMean(PAInst.Criteria == 0 & PAInst.Day == i  & RatSelection), 1);
end

for i=1:max(PAInst.Day(PAInst.Criteria == 0))
    PAInst.AverageCSMinCPLatMean(i) = nanmean(PAInst.CSMinCPLatMean(PAInst.Criteria == 0 & PAInst.Day == i & RatSelection));
    PAInst.Error_AverageCSMinCPLatMean(i) = nanste(PAInst.CSMinCPLatMean(PAInst.Criteria == 0 & PAInst.Day == i & RatSelection), 1);
end


% To find averages of CS+ and CS- Ratios by day
for i=1:max(PAInst.Day(PAInst.Criteria == 0))
    PAInst.AverageCSPlusRatio(i) = nanmean(PAInst.CSPlusRatio(PAInst.Criteria == 0 & PAInst.Day == i & RatSelection));
    PAInst.Error_AverageCSPlusRatio(i) = nanste(PAInst.CSPlusRatio(PAInst.Criteria == 0 & PAInst.Day == i & RatSelection), 1);
end

for i=1:max(PAInst.Day(PAInst.Criteria == 0))
    PAInst.AverageCSMinRatio(i) = nanmean(PAInst.CSMinRatio(PAInst.Criteria == 0 & PAInst.Day == i & RatSelection));
    PAInst.Error_AverageCSMinRatio(i) = nanste(PAInst.CSMinRatio(PAInst.Criteria == 0 & PAInst.Day == i & RatSelection), 1);
end

%% Save file
filename = 'InstTrainingAnalysis.mat';
save(filename);