如何從時間點爲numpy數組創建掩碼？

Question

數據是一個包含2500次測量時間序列的矩陣。我需要隨時間對每個時間序列進行平均，丟棄峯值附近記錄的數據點（在tspike-dt * 10 ... tspike + 10 * dt的時間間隔內）。針對每個神經元的尖峯時間數是可變的，並存儲在具有2500個條目的字典中。我當前的代碼遍歷神經元和尖峯時間，並將屏蔽值設置爲NaN。然後bottleneck.nanmean（）被調用。然而，這個代碼在當前版本中會變慢，我想知道更快的解決方案。謝謝！

import bottleneck 
import numpy as np 
from numpy.random import rand, randint 

t = 1 
dt = 1e-4 
N = 2500 
dtbin = 10*dt 

data = np.float32(ones((N, t/dt))) 
times = np.arange(0,t,dt) 
spiketimes = dict.fromkeys(np.arange(N)) 
for key in spiketimes: 
    spiketimes[key] = rand(randint(100)) 

means = np.empty(N) 

for i in range(N):   
    spike_times = spiketimes[i] 
    datarow = data[i] 
    if len(spike_times) > 0: 
    for spike_time in spike_times:       
     start=max(spike_time-dtbin,0) 
     end=min(spike_time+dtbin,t) 
     idx = np.all([times>=start,times<=end],0) 
     datarow[idx] = np.NaN 
    means[i] = bottleneck.nanmean(datarow) 

Answer 1

絕大多數的在你的代碼的處理時間來自這條線：

idx = np.all([times>=start,times<=end],0) 

這是因爲每個秒殺，你是在對抗開始和結束時間比較每個值。既然你有統一的時間步驟，在這個例子中（我想這是真的在你的數據也一樣），它是要快得多簡單計算的起始和結束的索引：

# This replaces the last loop in your example: 
for i in range(N):   
    spike_times = spiketimes[i] 
    datarow = data[i] 
    if len(spike_times) > 0: 
     for spike_time in spike_times: 
      start=max(spike_time-dtbin,0) 
      end=min(spike_time+dtbin,t) 
      #idx = np.all([times>=start,times<=end],0) 
      #datarow[idx] = np.NaN 
      datarow[int(start/dt):int(end/dt)] = np.NaN 
    ## replaced this with equivalent for testing 
    means[i] = datarow[~np.isnan(datarow)].mean() 

這減少了運行時間對我來說從大約100秒到大約1.5秒。 您還可以通過將spike_times上的循環向量化來削減更多時間。這樣做的效果將取決於您的數據的特徵（應該對高峯值速率最有效）：

kernel = np.ones(20, dtype=bool) 
for i in range(N):   
    spike_times = spiketimes[i] 
    datarow = data[i] 
    mask = np.zeros(len(datarow), dtype=bool) 
    indexes = (spike_times/dt).astype(int) 
    mask[indexes] = True 
    mask = np.convolve(mask, kernel)[10:-9] 

    means[i] = datarow[~mask].mean() 

Answer 2

而不是使用nanmean你可以只索引你需要的值，並使用mean的。

means[i] = data[ (times<start) | (times>end) ].mean() 

如果我誤解，你需要你的索引，你可以嘗試

means[i] = data[numpy.logical_not(np.all([times>=start,times<=end],0))].mean() 

你可能想不使用 if len(spike_times) > 0代碼

也（我假設你刪除尖峯時間在每次迭代或否則該語句將始終爲真，並且您將有一個無限循環），只能使用for spike_time in spike_times。