PySpike: 用于分析尖峰序列同步的开源Python软件包

⃝⃝可在www.sciencedirect.com上在线获取ScienceDirectSoftwareX 5（2016）183原始软件出版物www.elsevier.com/locate/softxPySpike-用于分析尖峰序列同步MarioMulansky，Thomas Kreuz复杂系统研究所，CNR，Via Madonna del Piano 10接收日期：2016年3月10日;接收日期：2016年7月14日;接受日期：2016年7月15日摘要了解大脑的功能是我们这个时代最大的挑战之一。实验记录的神经放电模式（尖峰序列）的分析在解决这个问题中起着至关重要的作用。这里介绍PySpike库，这是一个用于尖峰序列分析的Python包，提供了无参数和时间尺度独立的尖峰序列同步测量它允许计算相似性和相异性配置文件，平均值和距离矩阵。虽然PySpike主要专注于神经科学，但它也可以应用于气候研究或社会科学等其他领域。该软件包在Github和PyPI上以开源形式提供。c2016作者。由Elsevier B. V.发布。这是CC BY-NC-ND许可下的开放获取文章（http：//creativecommons. 组织/licenses/by-nc-nd/4. 0/）。关键词：同步; Spike train分析; Spike train距离; Python代码元数据当前代码版本v0.5.1用于此代码版本的代码/存储库的永久链接https://github.com/ElsevierSoftwareX/SOFTX-D-16-00032法律代码许可BSD许可使用git的代码版本控制系统使用Python、Cython的软件代码语言、工具和服务编译要求，操作环境依赖numpy，cython，matplotlib，nosetests如果可用，链接到开发人员文档/手册http://www.pyspike.de问题支持电子邮件mario.gmx.net1. 介绍深入了解大脑的内部运作仍然是一个很大程度上尚未解决的挑战，需要生物物理学，医学，实验以及计算神经科学的共同努力[1]。这一领域的科学进步的基础是神经活动的实验记录，通常以尖峰序列表示，即每个记录的神经元的尖峰时间列表。利用先进的现代记录技术，现在可以对神经活动进行高度并行的测量，通常会产生非常大的尖峰序列[2，3]。这将增加*通讯作者。电子邮件地址：mario. isc.cnr.it（M. Mulansky），thomas.cnr.it（T. Kreuz）。需要强大和高质量的数据分析工具，能够处理并行记录产生的大型数据集存在许多方法来分析尖峰序列数据，例如基于尖峰计数分布、尖峰间隔或精确的尖峰时间。一个非常重要的方法是量化尖峰序列之间的同步性。在过去的几十年中，已经提出了几种同步性测量[4-这里介绍的PySpike库1（图中所示的徽标）。 1）是一个Python包，允许一个人计算两个1www.pyspike.de。http://dx.doi.org/10.1016/j.softx.2016.07.0062352-7110/c2016作者。由爱思唯尔公司出版 这是一篇基于CC BY-NC-ND许可证的开放获取文章（http：//creativecommons. org/licenses/by-nc-nd/4.0/）。184M. Mulansky，T.Kreuz/SoftwareX 5（2016）183=={个联系我们Fig. 1. PySpike库的徽标。不同的相异度度量，ISI距离[9]，SPIKE距离[10]以及额外的相似性度量SPIKE同步[11，12]。这三种方法都是时间分辨的，无参数和时间尺度独立的，因此非常通用。例如，时间分辨意味着这些测量可以检测同步随时间的变化，而无参数使它们易于应用，具有明确的结果，因为不需要参数优化。这些措施已经在过去的许多实验研究中应用，例如[13虽然是在神经科学的背景下发展起来的，但这里讨论的同步性测量可以应用于由任何类型的事件序列组成的任何形式的离散时间序列。事实上，这些措施已经被用于其他几个研究领域，如气候研究[17]或社会科学[18，19]。PySpike是一个旨在使用Python脚本执行自动化数据分析的库。因此，它是SPIKY软件包的补充方法，2是一种用于棘波序列分析的Matlab框架，提供类似的功能，但还提供复杂的GUI [11，20]。最近开发了几个用于尖峰序列分析的其他软件包，特别是SyncPy3[21]，一个基于Python的GUI，用于量化时间序列中的同步性。但是，它目前不包括PySpike中实现的同步措施。在[22]中提出了尖峰序列距离测量4的C++实现，但它基于采样数据，因此性能相当差[11]。最后，作为[23]的一部分，还提供了一系列用于尖峰序列分析的科学软件，专门针对多变量记录。52. Spike Train距离这里，离散时间序列由尖峰序列表示，尖峰序列是表示事件（尖峰）在这些时间点的发生的时间点序列：时间分辨距离测量将一对尖峰序列s1、s2映射到轮廓s1、s2S（ t）上，其中0S（ t） 1。总距离值可以通过积分容易地获得：DS S（ t）dt。PySpike提供了三种这样的距离测量：ISI距离，SPIKE-distance和SPIKE-Synchronization。这些方法对锋电位序列的不同方面敏感2 http://www.fi.isc.cnr.it/users/thomas.kreuz/Source-Code/SPIKY.html网站。3 https://github.com/syncpy/SyncPy网站。4 https://github.com/modulus-metric网站。5 http://spiketrain-analysis.org/software网站。同步性（分别为尖峰间期、精确尖峰定时、尖峰匹配）。因此，方法的选择应该根据对信息如何在尖峰序列中编码的假设来确定。下面，我们简单介绍一下PySpike库中提供的措施。有关方法及其属性的详细讨论，请参见附录A和[12]。在[9]中引入的ISI距离分布I（ t）量化了两个锋电位序列的同时发放间期的相对差异。它本质上是测量两个锋电位序列瞬时速率的相对差异，但对精确的锋电位时间不敏感。ISI距离曲线是一个双变量分段常数函数。SPIKE-距离曲线S（ t），首次在[24]中引入并在[10]中细化，表示基于精确尖峰定时的相异度曲线。因此，尖峰距离根据与两个尖峰序列中的尖峰的精确重合的偏差来量化尖峰序列相异性。这导致双变量分段线性分布。虽然ISI和SPIKE距离曲线的基本定义是双变量（两个尖峰序列的距离曲线），但对多变量背景的概括是所有尖峰序列对的直接平均值[25]。SPIKE-Synchronization[11，12]是一种具有自适应符合窗口的直接、归一化符合计数器。它以两个尖峰序列之间的重合部分来量化相似性，因此是一个非常直观的度量。可以基于所有尖峰序列对来定义许多尖峰序列的尖峰同步的泛化，从而产生一致的多变量框架，其中相似性再次被量化为所有尖峰序列中的总体重合分数[11，12]。3. 封装结构3.1. 该尖峰序列PySpike库的中央数据结构是SpikeTrain，这是一个表示单个spike train的Python类这个类包含（排序的）尖峰时间作为一个numpy.array，以及尖峰序列的开始和结束时间。 这样的SpikeTrain对象可以通过提供尖峰时间直接创建，使用生成泊松尖峰函数从泊松过程随机生成，通过从txt加载尖峰序列从文本文件导入，或者通过从时间序列导入尖峰序列从时间序列导入。然后，这些对象用作计算距离测量的输入，参见。 图 二、3.2. 计算配置文件时间分辨是这里讨论的三个尖峰序列同步措施的主要优点。因此，PySpike包含计算同步曲线的功能：isi曲线计算分段常数ISI曲线I（t），尖峰曲线返回分段线性SPIKE曲线S（t），并且尖峰同步曲线产生离散的M. Mulansky，T.Kreuz/SoftwareX 5（2016）183185图二、PySpike包的结构SpikeTrain是中心类，提供了几个函数来计算相异性或相似性配置文件，以及距离/相似性值和矩阵。清单1：使用profile和distance函数。#g e n e r a r b i t r a r a r t r a i n s # w i t h s t a r t / end ti m e of 0和4st1= SpikeTrain（[1. 0，2。0，3。0]，edges=[0，4]）st2= SpikeTrain（[0. 五、三。0，3。5]，edges=[0，4]）st3= SpikeTrain（[2.五、三。8]，edges=[0，4]）#我爱你 p r o f i l e o f s p i k e T R A I N S1、2：isi_prof = isi_profile（st1， st2）#d i s t a n c e 作为一个v e r a g e o v er p r o f i l e伊西区=isi_prof. avrg（）# 如有疑问，请联系：直接控制从T H Es p i k e t r a i n sisi_dist = isi_distance（st1， st2）#p r o f i l e外来资产 3我想知道的是：spike_prof = spike_profile（[ st1， st2，st3 ]）#com p u t i n g p a i r w is s p i k e s y n c ma t i x #of the whole l i s to f s p i k e t a i ns：spike_trains = [ st1， st2， st3]M=spike_sync_matrix（spike_trains）# M I s一 3 x 3 numpy 阿雷阿雷SPIKE-同步配置文件Ck。 由于这三个配置文件具有不同的数学属性，因此它们在PySpike中由不同的Python类 表 示 ： PieceWiseConstFunc ， PieceWiseLinFunc 和DiscreteFunc，如图2所示。通过为每个对象引入特定的数据结构，PySpike可以利用它们的数学属性并使用高效的实现，同时仍然提供方便的用户界面。所有三个函数类都提供了一个avrg成员函数，该函数返回相应度量的在清单1中，计算并平均了几个此外，概要文件还提供了一个getplocomputdata成员函数，以方便概要文件的可视化，如清单2所示。3.3. 计算距离和距离矩阵除了计算时间分辨曲线，PySpike以及总的SPIKE-Synchronization值。相应的函数是isi距离、尖峰距离和尖峰同步，有关如何计算曲线的示例，请参见清单1，有关曲线图的示例，请参见图3。虽然也可以使用avrg（）成员函数从配置文件中计算整体距离，但特定的距离函数允许更有效的实现，因此具有更好的性能，如将在第4.2节中讨论的。当分析尖峰序列的集合时，人们可能不仅对整个集合的多变量距离感兴趣，而且对两个集合之间的成对距离感兴趣。所有的穗火车。为此，PySpike提供了isi距离矩阵，尖峰距离矩阵和尖峰同步矩阵函数。清单1显示了如何计算SPIKE-Synchronization矩阵的示例。由于PySpike中的所有距离度量都是对称的，因此距离矩阵也总是对称的。最后，请注意，距离和距离矩阵函数也可以用于计算选择性平均值，即分布图的指定间隔上的平均值 这是通过向距离函数提供（可选的）命名参数interval来实现的，对于配置文件的avrg成员函数也是如此。清单2给出了计算和绘制这种选择性距离矩阵的示例。这允许对轮廓进行非常详细的分析，例如通过比较重复刺激或触发之前和之后的平均距离图4示出了由具有不同平均间隔的距离矩阵揭示的变化的同步模式4. 实现和性能4.1. 代码库PySpike是按照现代编码标准和科学计算的最佳实践编写的，特别是遵循官方Python代码风格指南（PEP 8）。6可以还提供了直接计算总体ISI距离和SPIKE距离的距离值6https://www.python.org/dev/peps/pep-0008/。186M. Mulansky，T.Kreuz/SoftwareX 5（2016）183===清单2：加载加标序列并绘制曲线和距离矩阵。# l o a d s p i k e t r a i n D a t a2019 - 05 - 2501：0001：00“spike_trains.txt”，edges =（0，4000））# compute and plot I S I prof i l eisi_prof=isi_profile（spike_trains）x，y=isi_prof. get_plotmap_data（）plot（x，y，# p l o t 尖峰 D i s t a n c e 我是一个#w i t h s e l e c t i v e a v e r a g i n g t = 0.. . 1 0 0 0spike_mat 2019 - 05 - 22 00：01：00：00spike_trains，interval =[0，1000]）imshow（spike_mat）图三. M50峰电位序列的多变量ISI、SPIKE和SPIKE同步曲线（如上图所在该示例中，使用人工生成的尖峰序列，其中前半部分由叠加有具有越来越多抖动的几个同步事件的噪声尖峰组成后半部分包含越来越同步的事件，没有嘈杂的背景。来源：改编自[12]。与Python2和Python3运行时兼容，并作为BSD下的开源软件见图4。 一组M的尖峰距离矩阵40个人工道钉列车(top面板）具有变化的同步发射。锋电位序列分为四组，在观察间隔的前半部分和后半部分以两种不同的同步模式放电。通过对SPIKE距离矩阵（左侧T=[0，50]，右侧T=[50， 100]）进行选择性平均，可以清楚地捕捉到PySpike库在Github7上开发，并在PyPI上定期发布。8它使用nosetest框架进行了完全的单元测试，9包括持续集成服务器上的广泛测试矩阵，10使用Python版本2.6，2.7，3.3，3.4，3.5以及Python和Cython后端对每个提交执行所有测试。4.2. 性能随着神经科学实验能力的提高，现在可以在培养的神经元[30]和体内记录[2]中获得数千个和更多的尖峰序列的平行记录。这需要有效地实现尖峰序列分析工具，使人们能够在可接受的时间内处理如此然而，与C/C++等低级语言相比，Python的性能非常差因此，我们设计了PySpike库由两部分组成：（1）完全用Python实现的前端，代表PySpike功能与库用户的接口;（2）提供距离度量的数值实现的后端。对于后端，PySpike提供了两个版本，一个纯Python实现和一个基于Cython的更快的C实现。Cython后端是默认选择，只有当Cython不可用时，PySpike才回退到Python后端。利用Cython（以及C的速度），我们执照前端完全用Python实现，而对于关键的计算，Cython[27]版本提供了最大的性能，如下面的4.2作为大多数科学Python库，PySpike需要NumPy，因为numpy数组[28]是底层数据结构。此外，它被设计为与matplotlib交互，它可以绘制轮廓和距离矩阵，也可以很容易地与其他科学库（如scikit-learn）结合进行聚类分析[29]。能够获得超过一个因素的性能改进为200 μ图5显示了多变量分析的运行时的M1000泊松尖峰列车，每个包含N500尖峰，代表一个现实的实验情况。7 https://github.com/mariomulansky/PySpike/网站。8 https://pypi.python.org/pypi/pyspike/网站。9 https://nose.readthedocs.org/en/latest网站。10 https://travis-ci.org/mariomulansky/PySpike网站。M. Mulansky，T.Kreuz/SoftwareX 5（2016）183187==我∼I（t）= |−|我P我P我图五.用于计算一组平均包含N500个尖峰的M1000个尖峰序列的多变量ISI距离、尖峰距离和尖峰同步。深灰色条表示首先计算轮廓然后求平均的运行时，而浅灰色条表示距离值的直接计算。左侧面板显示Cython后端的运行时（以秒为单位），而右侧面板显示Python后端的运行时。对于两个后端，Cython（左）和Python（右），我们展示了ISI距离，SPIKE距离和SPIKE同步的完整多变量计算的运行时，首先计算轮廓，然后求平均值，以及直接距离函数。所有性能测量均在Intel Core i5- 3210 M CPU@2.50GHz上完成。首先，请注意，使用Python后端，计算时间达到几个小时，这显然超出了数据分析过程的任何可容忍的运行时间。然而，使用Cython后端，计算只需要10其次，如前所述，直接计算整体距离值（图5中的浅灰色）比首先计算轮廓然后求平均值（深灰色）要快得多，这就是PySpike提供此功能的原因。虽然图 5只显示了一种数据集大小的结果，即M=1000个尖峰序列，每个尖峰序列平均包含性能优异，因此适用于大型数据集的分析。除了定期维护和错误修复之外，PySpike未来的工作还将包括添加专门为突发尖峰序列设计的新同步措施[31];支持事件检测（从连续到离散数据的转换）;以及扩展尖峰序列生成功能。目前的工作重点是开发一个尖峰顺序指示器，允许通过量化尖峰序列集内的时间领导者/追随者属性来分析传播模式[32]。初步的Python实现已经存在，并将在不久的将来成为官方PySpike版本的一部分总之，我们相信PySpike库是一个有用的软件包，用于尖峰序列分析。它提供一个干净的，有良好记录的和一致的Python接口，用于计算尖峰序列同步的测量，同时提供出色的性能。致谢我们感谢N。 Bozanic和E. Ra？ isa？ nenformany有益的discussions. 此外 ，我 们感 谢我 。格 纳坚 科岛 Samuel 和R.Tomsett对PySpike的贡献M. M. 和T. K.感谢欧盟委员会通过玛丽·居里初始培训网络“神经工程转型技术（NETT）”项目289146提供的资金支持，以及T。K.玛丽·居里欧洲联合博士学位“复杂振荡系统：建模与分析（COSMOS）”项目642563也得到了支持此外，T.K.感谢意大利外交部对意大利-以色列神经科学联合实验室的支持。阑尾 数学定义A.1. ISI距离其中{t（1）}是第一尖峰序列的尖峰时间，其尖峰间间隔被给出为ν（1）=t（1）-t（1），并且类似地N超过500个尖峰，已知运行时线性i i+1i与尖峰的数量N成二次方关系，与尖峰序列的数量M成二次方关系（当对数增长M2时）。因此，较大（或较小）数据集的预期运行时间可以很容易地从图中的结果估计出来。 五、5. 结论和前景我们介绍了PySpike，这是一个Python库，用于测量实验和模拟尖峰序列数据的同步性。PySpike提供了三个无参数，时间尺度独立和多变量同步测量以及绘图，泊松尖峰序列生成和选择性平均的设施。它是在一个干净的、完全文档化的代码库对于第二道钉列（参见，图A.6）[9]。然后将ISI曲线计算为这些峰间期的归一化绝对差：ν（1）（t）ν（2）（t）max{ν（1）（t），ν（2）（t）}，t∈ [T0，T1]，（A.1）其中T0≤t（1，2）≤T1是尖峰序列的边缘A.2. SPIKE距离S（ t）的计算基于当前时间t周围的四个角尖峰：t（1）、（ 2）（t）和每个尖峰序列的后续尖峰t（1）、（ 2）（tP F遵循官方Python代码风格指南。此外，PySpike一体化 凭借Cython后端，它还显示（参见 图 [10]。对于每个角尖峰，计算到另一尖峰序列的最近尖峰的距离，例如：int（1）=min {|t（1）−t（2）|}，（A.2）188M. Mulansky，T.Kreuz/SoftwareX 5（2016）183- -2我=-我C（1）=MJ我JIj我我K我我Mk=1Mi j2i i−1jj−1其中，M是汇集的尖峰序列tk′中的尖峰总数，C表示重合尖 峰 的 总 数 与 其 他 两 种 测 量 方 法 相 比 ， SPIKE-Synchronization量化的是相似性而不是距离。引用图A.6. 计算ISI和SPIKE曲线所需的峰间期和时间差的局部定义而对于Lt（2）和Lt（1），（ 2）也是如此。这些距离是[1] QuianQuirogaR，Panzeri S. 神经编码原理 Boca Raton，FL，USA：CRC Taylor and Francis;2013.[2] Alivisatos AP等人，Nanotools for neuroscience and brain activitymapping。ACS Nano2013;7（3）：1850.[3] 大理石啊。等人，Physical Principles for Scalable Neural Recording，（137）第一百三十七章.[4] VictorJD，Purpura KP. 视觉皮层时间编码的性质和精确性：度量空间分析。神经生理学杂志1996;76：1310.[5] vanRossumMCW. 一种新颖的扣球距离。神经计算2001;13：751.[6] Quian Quiroga R，Kreuz T，Grassberger P.事件同步：测量同步性和时间延迟模式的简单快速方法。P F根据当前时间点t通过线性插值组合，得到以下量：t（1）（t） x（ 1）（t）+物理学修订版E 2002;66：041904。[7] Tiesing a PHE，Fellous JM，Sejnowski TJ.视皮层回路尖峰时间的调节。Nat Rev Neurosci2008;9：97.[8] 歌手W. 神经元中的分布式处理和时间编码S1（t）=PF Fν（1）（t）P，（A.3）netw orks. Cogn Neurodyn 2009;3：189.[9] Kreuz T、 莫尔曼 F、 安杰亚克 RG， 克拉斯科夫 一、LehnertzK，在耦合模型系统中测量同步：A并且类似地，S2（t）被定义用于第二尖峰序列。基于在这些量上，我们最终使用以下归一化得到SPIKE轮廓的定义：S1（t）ν（2）（ t）+S2（ t）ν（1）（ t）比较不同的方法。物理学D 2007;225：29。[10] [10]杨文辉，李文辉，李文辉. 监控尖峰序列同步。神经生理学杂志2013;109：1457。[11] [10] KreuzT，Mulansky M，Bozanic N. SPIKY：一个图形用户界面，用于监控尖峰序列同步。神经生理学杂志2015;113：3432。S（ t）=1（ν（1）（t）+ν（2）（t））2， t∈ [0，T].（A.4）[12] Mulansk yM，Bozanic N，Sburlea A，Kreuz T. 时间分辨指南和无参数测量的尖峰序列同步。在：基于事件的控制，通信和注意，这是分段线性函数， S1、2 是分段线性，而其他项是分段常数。A.3. SPIKE同步对于尖峰同步，为两个尖峰序列s（1）、（2）的每个尖峰定义重合指示符C（1）、（2），其中C如果t i处的尖峰是巧合的一部分，则为1，并且C i0如果没有[11]。巧合是根据适应性来定义的。符合窗口τ根据局部发射率：τ（ 1 ，2 ） =1min{ν（ 1 ） ， ν（ 1 ） ， ν（ 2 ） ， ν（ 2 ） } ，（ A.5）其中ν（1）、（ 2）是如上所述的峰间期。然后，重合指示符的值由下式给出：01- 02 - 02（|t（1）− t（2）|）<τ（1，2）信号处理（EBCCSP）。IEEE; 2015. p. 1.一、[13] Dodla R，Wilson CJ.耦合起搏神经元的非同步反应。物理学评论快报2009;102：068102。[14] Ibarz JM，Foffani G，Cid E，Inostroza M，de la Prida LM.癫痫海马快速波动的出现动力学。J Neurosci 2010; 30：16249.[15] 放大图片作者：WildieM. 建立神经元之间的通信通过竞争性夹带。Front Comput Neurosci 2012;5：62.[16] DiPoppaM，Gutkin BS. 背景活动控制中的相关性持久的状态稳定性，并允许执行工作记忆任务。 Front ComputNeurosci2013;7：139.[17] Malik N，Marwan N，Kurths J.南亚季风降水的空间结构和方向性。非线性过程地球物理2010;17：371.[18] W illiamsMJ，Whitaker RM，Allen SM. 测量个体规律性在人类拜访模式中。在：ASE社会计算国际会议论文集; 2012年。p.117.[19] Rabino witch TC，Knafo-Noam A.同步节奏的互动增强了儿童对彼此的相似性和亲密感。PLoS ONE2015;10：e0120878.0，否则。[20] [10]杨文，杨文.很机灵。Scholarpedia 2014;9（12）：32344. 版本号152246。以相同的方式，计算第二尖峰序列C（2）的重合指示符SPIKE-Synchronization配置文件然后由以下给出从合并的符合指标定义的crete函数{C k}={C（1）}<${C（2）}和尖峰时间：{t′}={t（1）}<${t（2）}。[21] [10]杨文，王文，王文. SyncPy：一个统一的开源同步分析库。在：在建模人际同步和影响力的第一次研讨会的会议记录。ACM; 2015. p. 41岁[22] Rusu CV，Florian RV.一类新的尖峰脉冲序列度量。神经计算2014;26（2）：306.[23] Gr unS，Rotter S. 并行spik e列车的分析Springer; 2010.[24] Kreuz T，Chicharro D，Greschner M，Andrzejak RG. 时 间 分辨由于积分不是离散函数定义的，过-所有的SPIKE同步值SYNC通过求和获得：SYNC=1Ck=C，（A.7）（A.6）M. Mulansky，T.Kreuz/SoftwareX 5（2016）183189和尖峰序列同步的时间尺度自适应措施。神经科学方法杂志2011;195：92.[25] KreuzT，Chicharro D，Andrzejak RG，Haas JS，Abarbanel HDI.测量多个尖峰序列的同步性。神经科学方法杂志2009;183：287.[26] Wilson G等，科学计算的最佳实践。PLoS Biol2014;12（1）：1。190M. Mulansky，T.Kreuz/SoftwareX 5（2016）183[27] [10]杨文军，李文军. Cython：两个世界中最好的。Comput SciEng2011;13（2）：31.[28] van der Walt S，Colbert SC，Varoquaux G. numpy数组：一种用于高效数值计算的结构。Comput Sci Eng2011;13（2）：22.[29] Pedre gosaF，et al. Scikit-learn：Python中的机器学习 J Mach LearnRes 2011;12：2825.[30] [10]张文辉，张文辉.有源像素传感器阵列，用于从单细胞到大规模神经元网络的高时空分辨率电生理记录。实验室芯片2009;9：2644.[31] [10]李晓，李晓.筹备中;2016年。[32] [10]李晓，李晓，李晓.筹备中;2016年。