Step-by-step calculation of rectangular window

조회 수: 6 (최근 30일)
Paul Hinze
Paul Hinze 2022년 7월 6일
댓글: Paul Hinze 2022년 7월 7일
Dear community!
I would like to make a step-by-step calculation of synaptic input. As you can see in the attached figure, the synaptic input rises only, but my goal is that after the given time constant tau, the synaptic input declines again, with the same value as it were raised.
The goal is:
A single spike should elicit a synaptic input that last for the given time of tau. Thus, the synaptic input should go down and not up.
Thank you in advance for your help.
%% Generate Single synaptic input shape
clear
close
clc
Ninput = 3;
Lrate = 20; % input intensity (lambda) [spikes/sec]
Tleng = 50; % simulated time length [ms]
dt = 0.01; % time step [ms]
% time vector
t = (0:round(Tleng/dt)-1)*dt; % time vector
tv = (0:round(Tleng/dt)-1)*dt; % time vector
% call poisson spike generator
spikes = poisson_train(Ninput, Lrate, Tleng, dt);
spikeall = sum(spikes,1); % total number of spikes at each time step
xalpha = 0;
yalpha = 0;
% creating output vector
g = zeros(1,length(spikeall));
% initial values
x = 0;
y = 0;
% factors used for calculation
Asyn = 3.5;
tau = 1; % ms
aa = Asyn; % amplitude factor for normalizing the input
ee = (dt/tau); % decrease factor used at each time step
% step-by-step calculation
for t = 1:length(spikeall)
y = y + aa*spikeall(t);
x = ee*x + dt*y;
g(t) = x; % store data
end
figure;
set(gcf, 'Position', [100 50 900 600]); % position and size
% plotting model synaptic input
subplot(2,1,1)
plot(tv,spikeall)
title('spikes')
subplot(2,1,2)
plot(tv,g);
title('simulated synaptic input');
xlabel('time [ms]');
  댓글 수: 2
Paul Hinze
Paul Hinze 2022년 7월 7일
At the end it should look like the second upper panel ;)
Paul Hinze
Paul Hinze 2022년 7월 7일
Function for spike generation:
function spikes = poisson_train(Ninput, Lrate, Tleng, dt)
% - Input:
% Ninput number of inputs (= number of source neurons or number of repetitions)
% Lrate input intensity (lambda) [spikes/sec]
% Tleng simulated time length [ms]
% dt time step [ms]
%
% - Output:
% spikes number of Poissonian spikes at each time step
% length(spikes) = round(Tleng/dt);
%
%% time vector
t = (0:round(Tleng/dt)-1)*dt; % time vector
%% parameters for Poisson process
ldt = Lrate*dt/1000.0; % (lambda) * (delta t) [no unit]
%% main part
% take random variable of length t
rd = rand(Ninput, length(t));
% compare the random variable with the probability to get spikes
spikes = ( rd < ldt ); % output spike vector

댓글을 달려면 로그인하십시오.

채택된 답변

Jonas
Jonas 2022년 7월 7일
i suggest you just add the rectangular window at the specified times given by the spices and the spice amplitude:
fs=10000;
t=0:1/fs:3;
tau=0.25; tauInSamp=tau*fs;
spikeTimes= (randi(t(end)*fs+1,10,1)-1); % generate 10 spikes somewhere in the time, unit is samples
spikeVals = randi(5,10,1);
resultingSpikes=zeros(size(t));
for spikeNr=1:numel(spikeTimes)
if spikeTimes(spikeNr)+tauInSamp-1<=numel(t) % if not near the end of the allowed time frame
resultingSpikes(spikeTimes(spikeNr):spikeTimes(spikeNr)+tauInSamp-1)=resultingSpikes(spikeTimes(spikeNr):spikeTimes(spikeNr)+tauInSamp-1)+spikeVals(spikeNr);
else
resultingSpikes(spikeTimes(spikeNr):end)=resultingSpikes(spikeTimes(spikeNr):end)+spikeVals(spikeNr);
end
end
ax1=subplot(2,1,1);
stem(spikeTimes/fs,spikeVals);
ax2=subplot(2,1,2);
plot(t,resultingSpikes);
linkaxes([ax1 ax2],'x')
  댓글 수: 2
Jonas
Jonas 2022년 7월 7일
note also the other, better readable solution further down the page
Paul Hinze
Paul Hinze 2022년 7월 7일
Thank you Jonas

댓글을 달려면 로그인하십시오.

추가 답변 (1개)

Jonas
Jonas 2022년 7월 7일
or you use the movsum function
fs=10000;
t=0:1/fs:3;
tau=0.25; tauInSamp=tau*fs;
spikeTimes= randi(numel(t),10,1); % generate 10 spikes somewhere in the time, unit is samples
spikeVals = randi(2,10,1)/2;
allSpikes=zeros(size(t));
allSpikes(spikeTimes)=spikeVals;
resultingSpikes=movsum(allSpikes,[tauInSamp-1 0]);
ax1=subplot(2,1,1);
stem(spikeTimes/fs,spikeVals);
ax2=subplot(2,1,2);
plot(t,resultingSpikes);
linkaxes()
  댓글 수: 1
Paul Hinze
Paul Hinze 2022년 7월 7일
Wow, even withut a loop!
This will speed up my code.
Thanks a lot ;)

댓글을 달려면 로그인하십시오.

카테고리

Help CenterFile Exchange에서 Just for fun에 대해 자세히 알아보기

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by