1

スパイク時間 (ニューロンからの活動電位) のベクトルと、刺激イベントのタイムスタンプのベクトルがあります。PSTH を作成して、刺激がニューロンのスパイク率に影響するかどうかを確認したいと考えています。各刺激イベントをループすることでこれを行うことができますが (以下の簡単な例を参照)、30,000 を超える刺激イベントがあり、多くのニューロンが記録されている長い実験では、これは非常に遅くなります。

forループなしでこれを行うにはどうすればよいですか?

遅い方法の例:

% set variables
spikeTimes = [0.9 1.1 1.2 2.5 2.8 3.1];
stimTimes = [1 2 3 4 5];        
preStimTime = 0.2;
postStimTime = 0.3;
for iStim = 1:length(stimTimes)
    % find spikes within time window
    inds = find((spikeTimes > (stimTimes(iStim) - preStimTime)) & (spikeTimes < (stimTimes(iStim) + postStimTime)));
    % align spike times to stimulus onset
    stimONtimes = spikeTimes(inds) - stimTimes(iStim);
    % store times in array for plotting
    PSTH_array(iStim,1:length(stimONtimes)) = stimONtimes;
end
4

2 に答える 2

0

最善の方法は、おそらく既存のヒストグラム関数を使用することです。それらは非常に高速で、必要なすべての情報を提供するはずです。もちろん、これはビンが重なっていないことを前提としています。サンプルデータを考えると:

spikeTimes = [0.9 1.1 1.2 2.5 2.8 3.1];
stimTimes = [1 2 3 4 5];        
preStimTime = 0.2;
postStimTime = 0.3;

次のようにビンを作成できます。

bins = sort([stimTimes - preStimTime, stimTimes + postStimTime])

また

bins = [stimTimes - preStimTime; stimTimes + postStimTime];
bins = bins(:).'

bins =
   0.80000   1.30000   1.80000   2.30000   2.80000   3.30000   3.80000   4.30000   4.80000   5.30000

次に、目的の結果と使用している MATLAB のバージョンに応じて、、、histcountsまたはdiscretizeを使用できます。histcを使用しますがhistc(それほど凝ったものは持っていないため)、入力は 3 つの関数すべてで同じです。histcounts(edges私たちには役に立たない)に対して 1 つの余分な出力があり、 discretize(実際のカウント) に対して 1 つ少ない出力があります。

[N, IDX] = histc(spikeTimes, bins)

N =    
   3   0   0   1   2   0   0   0   0   0

IDX =    
   1   1   1   4   5   5

ビンには(T(i) + postStimTime)との間の時間が含まれている(T(i+1) - preStimTime)ため、1 つおきのビンを取得する必要があります。

N = N(1:2:end)

N =
   3   0   2   0   0

同様に、奇数番号のタイムスロットで発生したスパイクにのみ関心があり、新しい に一致するようにインデックスを調整する必要がありますIDX

v = mod(IDX, 2)

v =
   1   1   1   0   1   1

IDX = ((IDX+1)/2).*v

IDX =
   1   1   1   0   3   3

結果は、最初に取得した結果と一致しています。ビン 1 に 3 つのスパイクがあり、ビン 3 に 2 つのスパイクがあります。

于 2016-08-14T22:44:26.670 に答える
0

以下は、すべてのスパイクと 2 つの仮定を 1 つのループで処理するソリューションです。

  • 刺激時間は一定間隔
  • 刺激間隔が PSTH 間隔より大きい

刺激時間が固定間隔であると仮定すると、次のようになります。

delta_times = mean(diff(stimTimes));
assert(max(abs(diff(stimTimes)-delta_times))<1e-3);

ここで、最初の刺激の前にスパイク時間を preStimTime に合わせます。

spikeTimes0 = spikeTimes - stimTimes(1) + preStimTime;

次に、2 番目の仮定を使用して、スパイクごとにカウントされる刺激を計算します。

assert((postStimTime-preStimTime)<dekta_times);
stimuli_index = floor(spikeTimes0 / delta_times); 

その刺激に対して相対的に計算します。

spike_time_from_stimuli = spikeTimes0 - stimuli_index*delta_times;

PSTH を 0.01 の精度で (他のすべての時間と同じ単位で) 構築しましょう。

dt = 0.01;
times_around_stimuli = preStimTime:dt:postStimTime;
n_time_bins = length(times_around_stimuli);
n_stimuli = length(stimTimes);
PSTH = zeros(n_stimuli, n_time_bins)
for i=1:length(spikeTimes)
     time_index = ceil(spike_time_from_stimuli(i) / dt);
     % Ignore time-bins far from the event
     if time_index > n_time_bins
           continue;
     end
     PSTH(stimuli_index(i),time_index) = PSTH(stimuli_index(i),time_index) + 1;
end
于 2016-08-14T20:11:55.153 に答える