如何从 numpy 数组的时间点创建掩码？

Question

data 是一个包含 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

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

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 .