matlab-(含教程)基于matlab的第三代SNN脉冲神经网络的仿真

clc; clear; close all; warning off; addpath(genpath(pwd)); %%%%%%%%%%%%%%%%%%%%% Set Paramaters and neuron model Hz = 1000; % 25ms ms = 1000/Hz; Ts = 1; T = 1000; % Enter T in miliseconds input = [0 0 1 1;... 0 1 0 1]; output = [0 0 0 1]; % The input neuron was a Poisson spike train at 40 Hz, % while the other were integrate-and-fire neurons. numNeuron = 60; numInput = numNeuron; NetArch = createNetwork(numInput, numNeuron); % Set paramaters of neurons Vth = 9; % Threshold voltage (mV) Reset = 0; % Reset Voltage (mV) RefPeriod = 2; % (ms) tau_mem = 2; C = 50; R = 1.2; NeuronModel = createNeuron(Vth, Reset, RefPeriod, tau_mem, C, R); clear numNeuron numInput Vth Reset RefPeriod tau_mem C R % hiddenNeuron = 40; inputNeuron = 2; Tr = .85; A = .5; % max signal for STDP tau = 100; % tau for STDP windowSize = 100; % for STDP %%%%%%%%%%%%%%%%%%%%% Create inputs and outpt spikes inputSpikes = zeros(2, T); outputSpikes = zeros(1, T); epoch = T/ms; for i = 1:floor(epoch/4) for j = 1:4 loc = fix(4*(i-1)*ms + find(input(1,:) == 1) * ms); inputSpikes(1, loc) = 1; loc = fix(4*(i-1)*ms + find(input(2,:) == 1) * ms); inputSpikes(2, loc) = 1; loc = fix(4*(i-1)*ms + find(output(1,:) == 1) * ms); outputSpikes(1, loc) = 1; end end figure, subplot(311), stem(inputSpikes(1, 1:400)); title('First input'); subplot(312), stem(inputSpikes(2, 1:400)); title('Second input'); subplot(313), stem(outputSpikes(1, 1:400)); title('Output'); %%%%%%%%%%%%%%%%%%%%% Create delays and wights % 40 neurons in input layer (20 neurons get spike from input 1, % 20 neurons get spike from input 2) % 40 neurons in hidden layer % 1 neuron in output layer delayInput = zeros(inputNeuron, hiddenNeuron); weightInput = rand(inputNeuron, hiddenNeuron); % delayHidden = zeros(hiddenNeuron, hiddenNeuron); % weightHidden = rand(hiddenNeuron, hiddenNeuron); delayOutput = zeros(hiddenNeuron, 1); weightOutput = rand(hiddenNeuron, 1); Itot = zeros(inputNeuron, T); for i = 2:T Itot(1, i) = (Itot(1, i-1)*Tr) + (NeuronModel.Vth*inputSpikes(1, i)); Itot(2, i) = (Itot(2, i-1)*Tr) + (NeuronModel.Vth*inputSpikes(2, i)); end figure, subplot(211), plot(Itot(1, 1:400)), title('First input'); subplot(212), plot(Itot(2, 1:400)), title('Second input'); % Input spike created and then with a delays return to hidden nerons %%%%%%%%%%%%%%%%%%%%% create hidden neuron for i = 1:hiddenNeuron NeuVar{i}.v = 0; NeuVar{i}.Spikes = zeros(1, T); NeuVar{i}.Refract = 0; %Refractory state of neuron NeuVar{i}.vStore = zeros(1, T); end %%%%%%%%%%%%%%%%%%%%% create output neuron NeuOut.v = 0; NeuOut.Spikes = zeros(1, T); NeuOut.Refract = 0; %Refractory state of neuron NeuOut.vStore = zeros(1, T); %%%%%%%%%%%%%%%%%%%%% Start periodic T ItotOut = zeros(hiddenNeuron, T); for t = 2:T disp(['time: ' num2str(t)]); for i = 1:hiddenNeuron % check if neuron fired in the previous time step if NeuVar{i}.v > NeuronModel.Vth NeuVar{i}.Spikes(t) = 1; NeuVar{i}.Refract = NeuronModel.RefPeriod; NeuVar{i}.v = NeuronModel.Reset; %if not then checks if neuron is in refractory period elseif NeuVar{i}.Refract > 0 NeuVar{i}.Refract = NeuVar{i}.Refract - Ts;%subtracts 1 time step from Refractory period NeuVar{i}.v = NeuronModel.Reset; %Sets mem voltage for this time step to 0 else %if neither of these then it calculates the new membrane voltage % signal from input neurons to hidden neurons I = Itot(:, t-1 - fix(delayInput(:,i))); x = weightInput(:, i) .* diag(I); deltav = sum(Ts .*( -1 .* NeuVar{i}.v / NeuronModel.tau_mem ... + NeuronModel.R_mem .* x / NeuronModel.tau_mem)); NeuVar{i}.v = NeuVar{i}.v + deltav; end NeuVar{i}.vStore(t) = NeuVar{i}.v; V(i, t) = NeuVar{i}.v; end for i = 1:hiddenNeuron ItotOut(i, t) = (ItotOut(i, t-1)*Tr) + (NeuronModel.Vth * NeuVar{i}.Spikes(t)); end % Output neuron if NeuOut.v > NeuronModel.Vth NeuOut.Spikes(t) = 1; NeuOut.Refract = NeuronModel.RefPeriod; NeuOut.v = NeuronModel.Reset; elseif NeuVar{i}.Refract > 0 NeuOut.Refract = NeuronModel.RefPeriod - Ts;%subtracts 1 time step from Refractory period NeuOut.v = NeuronModel.Reset; else I = ItotOut(:, t-1 - fix(delayOutput)); x = weightOutput .* diag(I); deltav = sum(Ts .*(((-1 .* NeuVar{i}.v / NeuronModel.tau_mem)) ... + ((NeuronModel.R_mem .* x)/NeuronModel.tau_mem))); NeuOut.v = NeuOut.v + deltav; end NeuOut.vStore(t) = NeuOut.v; % Update the output weigt out = outputSpikes(1:t); % Spikes of desired output in = NeuOut.Spikes(1:t); for i = 1:hiddenNeuron weightOutput(i) = STDP(in, out, A, A * 1.05, windowSize, tau, weightOutput(i)); end % Update the weight of input neurons for i = 1:hiddenNeuron % input 1 in = NeuVar{i}.Spikes(1:t); out = inputSpikes(1, 1:t); weightInput(1, i) = STDP(in, out, A, A * 1.05, windowSize, tau, weightInput(1, i)); % input 2 out = inputSpikes(2, 1:t); weightInput(2, i) = STDP(in, out, A, A * 1.05, windowSize, tau, weightInput(2, i)); end end figure, subplot(211), stem(NeuOut.Spikes), title('Spike of output input'); subplot(212), plot(NeuOut.vStore), title('V of output'); % Plot the spike that happened in hidden neurons for t = 2:T for i = 1:hiddenNeuron hiddenSpike(i, t) = NeuVar{i}.Spikes(t); end end plotMe(hiddenSpike)