Contents

clear;
close all;
clc;

Variable Initialization

Time variables

sweepNum = 1;
totalTime = 500;
deltaT = 0.05;
Rec_length = round(totalTime/deltaT);

% Size of Population
numNeurons = 200;

% Positional Info
coords = [];
xBounds = [-5,5];
yBounds = [0,4];

% Physiological Info
voltage = [];
spikes = [];
connection = zeros(numNeurons,numNeurons);
axons = [];
spikeMatrix = zeros(sweepNum, totalTime, numNeurons);
eTrain = zeros(Rec_length, numNeurons);
summedSpikeTrace = zeros(numNeurons,Rec_length); % will have the summed spike trains going into each neuron

% Analysis Info
corrs = [];

Generate Living Neurons

for i = 1:numNeurons
    % Generate position of neuron
    x = normrnd((xBounds(1)+xBounds(2))/2,1.5);
    y = normrnd((yBounds(1)+yBounds(2))/2,1);
    z = normrnd((yBounds(1)+yBounds(2))/2,1);
    coords = [coords; [x,y,z]];

    % Generate activity of neuron
    %V = HHspike(0.01);
    V = synBarrage(100, totalTime, sweepNum, deltaT);
    % Add to array of voltages
    voltage = [voltage; V];
    [~,spikeTrace] = countSpikes(V);
    spikes = [spikes; spikeTrace];
end

Generate Connections

for j = 1:numNeurons-1


    for i = 1:numNeurons-1
        % Generate Connections
        prob = normrnd(0.5,0.25);
        if (i==j) || (connection(i,j)==1)
            dist = 100000000; % don't let it synapse onto itself or let back-propogation
        end
        dist = pdist([coords(j,:); coords(i,:)],'Euclidean'); % get Euclidean distance
        if ((1/(2*dist))*prob > 0.5) && ((1/dist) ~= Inf) % if the generated probability times a distance weighting factor is over 0.5...
            connection(j,i) = 1; % ...there is a connection

            summedSpikeTrace(j,:) = summedSpikeTrace(j,:) + spikes(i,:); % add up the spikes
            % slower and unnecessary but easier for me to keep track of coordinates
            x1 = coords(j,1);
            y1 = coords(j,2);
            z1 = coords(j,3);
            x2 = coords(i,1);
            y2 = coords(i,2);
            z2 = coords(i,3);

            axons = [axons; [x1 x2],[y1 y2], [z1 z2]];
        end
    end
end

for i = 1:numNeurons
    alpha = genPSC(summedSpikeTrace(i,:)); % generate a post-synaptic current from the input spikes
    convolved = conv(summedSpikeTrace(i,:),alpha); % convolve the input spikes with generated alpha function
    convolved(Rec_length+1:Rec_length+Rec_length-1) = []; % clean up the edges
    voltage(i,:) = HHspike(convolved); % use the product of the convolution as input current to the receiving HH unit
end

Visualization

figure(1)
scatter3(coords(:,1),coords(:,2),coords(:,3),'ro');
hold on
scatter3(coords(:,1),coords(:,2),coords(:,3),'r+'); scatter3(coords(:,1),coords(:,2),coords(:,3),'rx');
for i = 1:size(axons,1)
    plot3([axons(i,1) axons(i,2)],[axons(i,3) axons(i,4)],[axons(i,5) axons(i,6)], 'r');
end

%figure(2)
%for x = 1:numNeurons
%    plot(voltage(x,:));
%    hold on
%end

%figure(3)
%hist(corrs');

Published with MATLAB® R2018a