%% Exploring individual gene expression
%% Load data

clearvars
clc
load filt_neurons-bc-somas-proj-axonal-clust_nocompr.mat

%


output_dir=fullfile('figs1','gene_expression_plots');
mkdir(output_dir)

%
%% Plot the expression patterns of individual genes
% Optional task: play with what genes to visualize.

genes_to_plot = [53 79 92]; % pick the genes you want to plot here

slices = unique(filt_neurons.slice);
logcp1000=log2(1+filt_neurons.expmat./sum(filt_neurons.expmat,2)*1000); % play with the scaling factor to change how the plot looks

for g = 1:numel(genes_to_plot)
    figure('Position',[70 70 1200 400]);

    for s= 1:numel(slices)
    subplot(1,numel(slices),s);
        in_slice=filt_neurons.slice==slices(s);

        %plot all neurons, colored by the logcp100 counts for the gene of
        %interest
        scatter(filt_neurons.pos(in_slice,1), ...
            filt_neurons.pos(in_slice,2), ...
            2, ...
            logcp1000(in_slice,g),...
            'filled' ...
            );
        xl=xlim;
        yl=ylim;
        set(gca, ...
            'ydir','reverse' ...
            );
        axis off
        pbaspect([range(xl),range(yl),1]);
        colorbar;
        colormap(sky)
        title(['Slice ',num2str(s)]);
    end
    sgtitle(['Gene ',filt_neurons.genes{genes_to_plot(g)}])
    exportgraphics(gcf,fullfile(output_dir,[filt_neurons.genes{g},'.pdf']),'ContentType','Image','Resolution',300);
end