{{tag>basal_ganglia computational_model decision_making}}

=====JC/hackathon for Gurney et al., 2001a/b=====

A computational model of action selection in the basal ganglia. I. A new
functional anatomy.

Gurney K, Prescott TJ, Redgrave P.

http://www.ncbi.nlm.nih.gov/pubmed/11417052

(see also the companion modeling paper,
http://www.ncbi.nlm.nih.gov/pubmed/11417053)


====Model information====
on ModelDB: [[http://senselab.med.yale.edu/modeldb/ShowModel.asp?model=83560]]

(This page contains a zip file with the Matlab scripts and some READMEs)

===Reproducing Gurney 2001a Fig3b===
Our goal with this script is to reproduce Fig3b:
{{ :jclub:fig3b.png?600x300 |}}

To recreate the input, activation and output values from Fig.3b we first need to set our parameters based on those found in the paper.

<code matlab>
  %% Parameters
  params.numChans = 4;                 % number of channels in the paper
  params.posW = 0.45;                  % The positive weights between input-output nodes
  params.negW = -1.35;                 % The negtive weights between input-output nodes
  params.chanIDs = 1:params.numChans;  % channel IDs (used for plots)
  W = ones(params.numChans);           %synaptic weight matrix (inputs-->outputs only)
  epsilon = -0.1;  % 
  m = 1;
</code>

Now we need to build the weight and input matrices

<code matlab>
%% weight matrix
weights = W*params.posW;
weights(eye(params.numChans)~=0)=params.negWeight; %find the diagonals and replaces them with the negative weights

%% input matrix
x = [.27, .77, .35, .83];    % from input fig3b
</code>

We can now take the input values and use a matrix multiplication to combine them with the weights to get the activation values for each channel

<code matlab>
%% activation
a = x*weights;
</code>

With the activation values for each channel we can then pass them through a piecewise linear squashing function:
<code matlab>
function [y] = piecewise(a, epsilon, m)
% piecewise
%
% Takes a matrix of size n a as the activation values
% a(a<epsilon) = 0;
% a((epsilon <= a) & (a <= (1/m+epsilon))) = m*(a-epsilon);
% a(a>(1/m+epsilon)) = 1;

for i = 1:length(a)
    if a(i) < epsilon;
        a(i) = 0;
    elseif ((epsilon <= a(i)) && (a(i) <= (1/m+epsilon)));
        a(i) = m*(a(i)-epsilon);
    else a(i) > 1/m+epsilon;
        a(i) = 1;
    end
    
end

y = a;  % corresponds to the output of the nodes 

end
</code>

We now have an input(x), activation(a) and output(y) values for each of the channels. Using the subplot functions we can put all three together to reproduce Fig3b.

<code matlab>
%% plot the results
figure(1)
%plot the input
subplot(1,3,1)                   % subplot(#columns, #rows, plot index)
bar(params.chanIDs, x ,0.5, 'k') % bar(x value, y value, width, color)
title('Input', 'FontSize', 24)
xlabel('Channel', 'FontSize', 24)
ylabel('Signal Level', 'FontSize', 24)
ylim([0 1])
set(gca,'YTick',-0:.25:1)      % sets the major ticks for a certain range

%plot the activation level
subplot(132)
bar(params.chanIDs, a,0.5, 'k')
title('Activation', 'FontSize', 24)
xlabel('Channel', 'FontSize', 24)
ylim([-.5 .5])
set(gca,'YTick',-0.5:.25:0.5)

% plot the output
subplot(133)
bar(params.chanIDs, y,0.5, 'k')
xlabel('Channel', 'FontSize', 24)
title('Output', 'FontSize', 24)
ylim([0 1])
hold on 
line(0:0.005:params.numChans+1, 0.1,'Color','k') % adds a line at yo like Fig3b
set(gca,'YTick',-0:.25:1)
text(0.1, 0.09, 'yo')               % adds text 'yo' near the line. 
</code>

This should yield a similar figure to Figb:
{{ :jclub:gurneyfig3b.png?600x300 |}}

Which looks a lot like Fig3b