n_regions=10;

show_simultaneous_results= false;

% prior distribution of location of rat
% it must be a vector of length n_regions
P_X = [4, 1, 1, 1, 1, 1, 1, 1, 1, 4];

% prior distribution must sum to 1, so normalize by its sum
P_X = P_X/sum(P_X);

% conditional distribution of neuron 1 (R1=1) given location X
% it must be a vector of length n_regions
P_R1_given_X = [0.05, 0.05, 0.1, 0.2, 0.15, 0.05, 0.05, 0.05, 0.05, 0.05]; 

% conditional distribution of neuron 2 (R2=1) given location X
% it must be a vector of length n_regions
P_R2_given_X = [0.01, 0.01, 0.01, 0.01, 0.04, 0.06, 0.03, 0.01, 0.01, 0.01];

% conditional distribution simultaneousl spikes (R1=1 and R2=1) given location X
P_R1_R2_given_X = P_R1_given_X.*P_R2_given_X;

% marginal distribution of spiking P(R1=1), P(R2=1), and P(R1=1, R2=1)
P_R1 = P_X * P_R1_given_X';  % dot product
P_R2 = P_X * P_R2_given_X';
P_R1_R2 = P_X * P_R1_R2_given_X';


% Bayes theorem for conditional distribution of location given spiking
P_X_given_R1 = P_R1_given_X.*P_X/P_R1;
P_X_given_R2 = P_R2_given_X.*P_X/P_R2;
P_X_given_R1_R2 = P_R1_R2_given_X.*P_X/P_R1_R2;


% calculate max probability for upper limit of y-axes
maxp = max([P_X, P_X_given_R1, P_X_given_R2, P_X_given_R1_R2])*1.05;

% calculate max spike probability for upper limit of y-axes
maxsp = max([P_R1_given_X, P_R2_given_X])*1.05;

% values of X for the x-axis
X=1:n_regions;

% colors in rgb
darkgreen = [0, 0.5, 0];
darkmagenta = [0.4, 0, 0.4];
darkblue = [0, 0, 0.5];

figure(1)
clf

subplot(3,2,1)
bar(X, P_R1_given_X, 'facecolor', darkgreen);
xlabel('X')
ylabel('Spiking probability')
title('Probability of neuron 1 spike')
axis([0, n_regions+0.9, 0, maxsp]);

subplot(3,2,2)
bar(X, P_R2_given_X, 'facecolor', darkgreen);
xlabel('X')
ylabel('Spiking probability')
title('Probability of neuron 2 spike')
axis([0, n_regions+0.9, 0, maxsp]);


% bar plot of prior probability distribution of location
subplot(3,2,3)
bar(X,P_X, 'facecolor', darkblue)
xlabel('X')
ylabel('Probability')
title('Prior distribution')
axis([0, n_regions+0.9, 0, maxp]);

if(show_simultaneous_results)
    % bar plot of probability distribution of location given 
    % 1 spike in both neurons measured
    subplot(3,2,4)
    bar(X, P_X_given_R1_R2, 'facecolor', darkmagenta)
    xlabel('X')
    ylabel('Probability')
    title('Post. dist. given R_1=1 and R_2=1')
    axis([0, n_regions+0.9, 0, maxp]);
end

% bar plot of probability distribution of location given
% a spike in neuron 1 was measured
subplot(3,2,5)
bar(X, P_X_given_R1, 'facecolor', darkmagenta)
xlabel('X')
ylabel('Probability')
title('Post. dist. given R_1=1')
axis([0, n_regions+0.9, 0, maxp]);

% bar plot of probability distribution of location given
% a spike in neuron 2 was measured
subplot(3,2,6)
bar(X, P_X_given_R2, 'facecolor', darkmagenta)
xlabel('X')
ylabel('Probability')
title('Post. dist. given R_2=1')
axis([0, n_regions+0.9, 0, maxp]);


