imbal_plots.m

function [eligible_people, high_accs] = imbal_plots(dirstring,...
    beh_type, plotting_param, filepath, file_ext, Xc, yc, ynan, p, pc, c)
  
    people = 1:136;
    include = zeros(max(size(people)),1);

    % Change the plot titles
    if contains(file_ext,'yes')
        kstring = "Direct Kernel";
    else
        kstring = "Linear";
    end
    
    FigH = figure('Position', get(0, 'Screensize'));
    person_i = 0;
    eligible_people = 1:136;
    highest_number = sum(p == mode(p(ismember(p,eligible_people)))) + 10;


    hex = ['#9ea2a2','#d6001c']';
    map = sscanf(hex','#%2x%2x%2x',[3,size(hex,1)]).' / 255;
    mC = [255 255 255];
    for	r=1:length(eligible_people)
        person                      =  eligible_people(r);
        index                       =  pc == person;
        Xp                          =  Xc(index,:);
        yp                          =  yc(index,1);

        index_1                     =  yp == 1;
        Xp1                         =  Xp(index_1,:);
        index_2                     =  yp == 0;
        Xp2                         =  Xp(index_2,:);
        Xp_0                        =  [Xp1;Xp2];
        Np_1                        =  size(Xp1,1);
        Np_2                        =  size(Xp2,1);
        yp_0                        =  [ones(Np_1,1);zeros(Np_2,1)];
        K                           =  [Np_1 Np_2];
        

        load(strcat(filepath,beh_type,file_ext,string(person),'.mat'),'WT1');



        if WT1{1}.Best == -1
            continue
        elseif sum(index) == 0
            continue
        else
            include(person) = 1;
        end

        BCRs = NaN(30,1);
        for i=1:30
            BCRs(i) =  WT1{i}.BCR;
        end

        if sum(isnan(BCRs)) > 10 % Ensure that there are enough BCRs in the model
            include(person) = 0;
            continue
        end
        


%         Start_List = cellfun(@(c) [c.Best], WT1, 'UniformOutput', false);
%         test = sum(cell2mat(Start_List) == 1);
    %     if test < 20
    %         error('test less than 30') % double check enough BCRs
    %     end


        person_i = person_i + 1;
        % Get the amount of class imbalance for each person
        all_imbals(person_i) = sum(index_1) / length(yp);
        if or(all_imbals(person_i) > 0.9, all_imbals(person_i) < 0.1)
            all_imbals(person_i)
            error("class imbal problem")
        end
        
        
        ydata = ynan(p==person)';
        ydata = [ydata NaN(1,highest_number - length(ydata))];
        hmap(person_i,:) = ydata;    
        outfile(:,person_i) = BCRs;
        
        BCRs = NaN(30,1);
        SENs = NaN(30,1);
        SPEs = NaN(30,1);
        P = 0;
        N = 0;
        avgs = zeros(30,1);
        % Iterate over all 30 random seeds

        for i=1:30
            currP = WT1{i}.TP + WT1{i}.FN;
            currN = WT1{i}.TN + WT1{i}.FP;
            test_set = currP + currN;
            curr_avg = currP / (currP + currN);
            avgs(i) = curr_avg;
            P = P + currP;
            N = N + currN;
            if or(test_set == currP, test_set == currN)
                continue % There is not a classification problem
            end
            BCRs(i) =  WT1{i}.BCR;
            SENs(i) =  WT1{i}.SEN;
            SPEs(i) =  WT1{i}.SPE;

        end


        tot = P + N;
        [mean(avgs), P/tot];
        sens_avg(person_i) = nanmean(SENs);
        bcrs_avg(person_i) = nanmean(BCRs);
        spes_avg(person_i) = nanmean(SPEs);
        
    end

    plot_out.SEN = sens_avg;
    plot_out.SPE = spes_avg;
    plot_out.BCR = bcrs_avg;
%     
%     disp(max(all_imbals(include==1)))
%     disp(min(all_imbals(include==1)))

    tiledlayout(1,1,'TileSpacing','Compact','Padding','Compact');
    nexttile;
    frac = 20;
    
    % Create a fraction for the main plot
    h = heatmap(hmap,'Colormap',map, ...
            'MissingDataColor', mC./255, ...
            'ColorbarVisible','off','FontSize',18);
    h.GridVisible = 'off';
    my_str = '';
    
    % Figure lettering
    if strcmp(beh_type,'AGG')
        my_str = 'a';
    elseif strcmp(beh_type, 'SIB')
        my_str = 'b';
    elseif strcmp(beh_type, 'BOTH')
        my_str = 'c';
    end
    
    % Label graph
    title(strcat("(", my_str, ") ", beh_type, " Behavior Patterns over Time"))
    xlabel('Days')
    ylabel('Person Number')
    s = struct(h);
    s.XAxis.TickLabelRotation = 45;
    XLabels = 1:highest_number;
    CustomXLabels = string(XLabels);
    CustomXLabels(mod(XLabels,50) ~= 0) = " ";
    h.XDisplayLabels = CustomXLabels;

    YLabels = 1:sum(include);
    CustomYLabels = string(YLabels);
    CustomYLabels(mod(YLabels,10) ~= 0) = " ";
    h.YDisplayLabels = CustomYLabels;
    fig = 0;
    set( findall(fig, '-property', 'fontname'), 'fontname', 'Palatino Linotype');
    saveas(FigH, strcat(dirstring,'/',beh_type,c,'_behavior_patterns.png'),'png');

    % Create a fraction of subplot for the legend
    FigH = figure('Position', get(0, 'Screensize'));
    q = subplot(frac,1,1);
    h = zeros(3, 1);
    h(1) = plot(NaN,NaN,'Color',map(2,:),'LineWidth',10);
    hold on;
    h(2) = plot(NaN,NaN,'Color',map(1,:),'LineWidth',10);
    h(3) = plot(NaN,NaN,'Color','white','LineWidth',10);
    set(q,'Visible','off');
    legend(h, 'Behavior','No Behavior','No Data','Location','northeastoutside','FontSize',14);
    saveas(FigH, strcat(dirstring,'/',beh_type,c,'_legend.png'),'png');
    
    
    FigH = figure('Position', get(0, 'Screensize'));
    tiledlayout(1,1,'TileSpacing','Compact','Padding','Compact');
    nexttile;
    ylim([0 1]);
    ax = gca; 
    ax.FontSize = 18;
    ax.YLim = [0 1];
    boxplot(ax, outfile,'Colors','k')
    ylim([0 1]);
    ylim([0 1]);
    xticks([0:5:80]);

    ax = gca; 
    ax.FontSize = 18;
    xticklabels(string([0:5:80]));
    %yticklabels(string([0:0.1:1]))
    xtickangle(45);

    % Title and label figure
    title(strcat(plotting_param, " for ",beh_type," Individuals, ", kstring," with Logistic Activation"), 'FontSize',20)
    fig = 0;
    set( findall(fig, '-property', 'fontname'), 'fontname', 'Palatino Linotype');
    ylabel(plotting_param, 'FontSize',20);
    xlabel('Person Number','FontSize',20);
    saveas(FigH, strcat(dirstring,'/',beh_type,c,'_kernel_boxplot.png'),'png');

    imbals = all_imbals;
    for i=1:length(imbals)
        max_imbal(i) = max([imbals(i) 1-imbals(i)]);
    end
    plotting_params = ["SEN","SPE","BCR"];
    tiledlayout(1,3,'TileSpacing','Compact','Padding','Compact');
    
    for j=1:length(plotting_params)
        nexttile;
        plotting_param = plotting_params(j);

        plot(imbals, plot_out.(plotting_param), 'r.','MarkerSize',15);
        ax = gca; 
        ax.FontSize = 14; 

        
        xlabel(strcat("% ", beh_type," Behavior"),'FontSize',14)
        title(strcat("% ",beh_type," Behavior vs. Avg. ", plotting_param),'FontSize',14)
        xlim([0 1]);
        ylim([0 1]);
        if strcmp(plotting_param, 'BCR')
            alpha(0.1)
            hold off;
            toplot = plot_out.(plotting_param);
            plot(max_imbal, toplot, 'r.','MarkerSize',15);
            alpha(0.1)
            hold on;
            plot(max_imbal(toplot > 0.8), toplot(toplot > 0.8), 'r.','MarkerSize',15,'LineWidth',1,'MarkerEdgeColor','black');
            alpha(0.1)
            xlim([0.5 1]);
            ylim([0.5 1]);
            hold on; 
            yline(0.8,'--','Color',[0.4 0.4 0.4]);
            xlabel('Max Class Imbal.','FontSize',14);
            title(strcat(beh_type," Class Imbal. vs. Avg. ", plotting_param),'FontSize',14);
        end
        hold on;
        if strcmp(plotting_param, 'SPE')
            plot(linspace(0,1,100),linspace(1,0,100),'k-');
        else
            plot(linspace(0,1,100),linspace(0,1,100),'k-');
        end
        ylabel(strcat("Avg. ",plotting_param),'FontSize',14);
        ax = gca; 
        ax.FontSize = 14; 
        axis square;

    end
    
    
    fig = 0;    
    set( findall(fig, '-property', 'fontname'), 'fontname', 'Palatino Linotype');
    saveas(FigH, strcat(dirstring,'/',beh_type,c,plotting_param,'_dotplot.png'),'png');
    eligible_people = find(include==1);
    save(strcat(dirstring, '/',beh_type,'people2'),"eligible_people");

    high_accs = bcrs_avg > 0.8;
    save(strcat(dirstring, '/',beh_type,'people2_highacc.mat'),"high_accs");
    fig = 0;   
    set( findall(fig, '-property', 'fontname'), 'fontname', 'Palatino Linotype');
end