CN106716848B

CN106716848B - 用于压缩测序数据的方法、系统和计算机可读媒体

Info

Publication number: CN106716848B
Application number: CN201580029528.4A
Authority: CN
Inventors: B.唐内特
Original assignee: Life Technologies Corp
Current assignee: Life Technologies Corp
Priority date: 2014-06-04
Filing date: 2015-06-03
Publication date: 2020-10-30
Anticipated expiration: 2035-06-03
Also published as: EP3152839B1; CN106716848A; US10254242B2; WO2015187832A1; EP3152839A4; US20150355138A1; EP3152839A1

Abstract

本发明公开用于压缩测序数据的方法、系统和计算机可读媒体。一种方法包含：接收与传感器阵列上发生的化学事件相关联的波形数据，所述波形数据包含所述传感器阵列的对应多个位置的多个基于时间的波形；通过至少一个处理器将所述波形数据的每一基于时间的波形转换成频域谱；通过所述至少一个处理器基于所述多个频域谱产生关键帧；通过所述至少一个处理器针对所述频域谱中的每一个计算所述频域谱与所述关键帧之间的差；以及通过所述至少一个处理器编码所述频域谱与所述关键帧之间的每一经过计算的差。

Description

用于压缩测序数据的方法、系统和计算机可读媒体

相关申请案的交叉参考

本申请案依据35 U.S.C. §119(e)主张2014年6月4日申请的第62/007,435号美国临时申请案的优先权，所述美国临时申请案的标题为“用于压缩测序数据的方法、系统和计算机可读媒体（Methods, Systems, and Computer-Readable Media for Compression ofSequencing Data）”，且前述申请案的内容以全文引用的方式并入本文中。

技术领域

本发明大体上是针对关于压缩测序数据的发明方法、系统和计算机可读媒体，所述测序数据是通过检测和测量一个或多个分析物而获得的，所述一个或多个分析物包含与核酸合成反应相关联或由核酸合成反应产生的分析物。

背景技术

电子装置和组件已在包含化学和生物学的生命科学中找到众多应用，尤其是用于检测和测量各种化学和生物反应以及识别、检测和测量各种化合物。一个此类电子装置被称作离子敏感场效应晶体管（“ISFET”）。ISFET促进对溶液的氢离子浓度（通常表示为“pH”）的测量。

更具体来说，ISFET为以类似于金属氧化物半导体场效应晶体管（“MOSFET”）的方式操作的阻抗变换装置，且确切地说，经配置以选择性地测量溶液中的离子活性（例如，溶液中的氢离子为“分析物”）。

随着取样数据速率变快以及ISFET的传感器阵列密度变高，可生成大量数据。因此期望降低存储器消耗，同时维持数据的质量。下文详细论述的至少某些方法的目标尤其在于准确地采集与生物/化学事件相关联的数据，同时降低与数据相关联的噪声。此目标可通过实施下文描述的压缩技术来达成。结果，可减少存储的数据量。

发明内容

实施例公开用于压缩测序数据的方法、系统和计算机可读媒体。

根据某些实施例，公开用于压缩测序数据的计算机实施的方法。一种方法包含：接收与传感器阵列上发生的化学事件相关联的波形数据，所述波形数据包含传感器阵列的对应多个位置的多个基于时间的波形；通过至少一个处理器将波形数据的每一基于时间的波形转换成频域谱；通过至少一个处理器基于多个频域谱产生关键帧；通过至少一个处理器针对频域谱中的每一个计算频域谱与关键帧之间的差；以及通过至少一个处理器编码频域谱与关键帧之间的每一经过计算的差。

根据某些实施例，公开用于压缩测序数据的系统。一个系统包含存储用于压缩测序数据的指令的数据存储装置；以及经配置以执行指令以执行方法的处理器，所述方法包含：接收与传感器阵列上发生的化学事件相关联的波形数据，所述波形数据包含传感器阵列的对应多个位置的多个基于时间的波形；将波形数据的每一基于时间的波形转换成频域谱；基于多个频域谱产生关键帧；针对频域谱中的每一个计算频域谱与关键帧之间的差；以及编码频域谱与关键帧之间的每一经过计算的差。

根据某些实施例，公开存储指令的非暂时性计算机可读媒体，所述指令在由计算机执行时使得计算机执行压缩测序数据的方法。一个计算机可读媒体包含如下方法：接收与传感器阵列上发生的化学事件相关联的波形数据，所述波形数据包含传感器阵列的对应多个位置的多个基于时间的波形；通过至少一个处理器将波形数据的每一基于时间的波形转换成频域谱；通过至少一个处理器基于多个频域谱产生关键帧；通过至少一个处理器针对频域谱中的每一个计算频域谱与关键帧之间的差；以及通过至少一个处理器编码频域谱与关键帧之间的每一经过计算的差。

所公开的实施例的额外目的和优点将在以下描述中部分阐述，且部分将从所述描述显而易见，或可通过实践所公开的实施例习得。所公开的实施例的目的和优点将借助于所附权利要求书中特别指出的元件和组合来实现和获得。

应理解，前文总体描述以及以下详细描述都仅仅是示范性以及说明性的，并且不限制所公开的实施例的范围，如权利要求书所阐述。

附图说明

并入在本说明书中并且构成本说明书的一部分的附图说明了各种示范性实施例，并且与所述描述一起用于解释所公开的实施例的原理。

图1描绘根据本发明的示范性实施例的包含大规模chemFET阵列的核酸处理系统；

图2描绘根据本发明的示范性实施例的chemFET传感器阵列的示范性CMOS IC芯片实施方案的框图；

图3描绘根据本发明的示范性实施例的时域波形；

图4描绘根据本发明的示范性实施例的频域波形频谱；

图5描绘根据本发明的示范性实施例的传感器阵列中的所有井的平均“均值”频谱；

图6描绘根据本发明的示范性实施例的图5中所示的平均“均值”频谱的截断平均“均值”频谱；

图7描绘根据本发明的示范性实施例的平均“均值”频谱（关键帧）与个别井之间的频谱相关；

图8描绘根据本发明的示范性实施例的基于频谱相关拒绝的井迹线；

图9A到9I描绘根据本发明的示范性实施例的平均“均值”频谱与个别井之间的各种示范性频谱相关；

图10描绘根据本发明的示范性实施例的归一化频率分量熵；

图11描绘根据本发明的示范性实施例的每一频率分量的位；

图12A和12B分别描绘根据本发明的示范性实施例的示范性关键帧和井的频谱；

图13描绘根据本发明的示范性实施例的来自关键帧的示范性井频谱增量；

图14描绘根据本发明的示范性实施例的压缩方法；

图15描绘根据本发明的示范性实施例的经重建构井频谱；

图16描绘根据本发明的示范性实施例的经重建构时域波形；

图17A到17J描绘根据本发明的实施例的传感器阵列的各种井的示范性原始时域波和经重建构时域波形；

图18A和18B描绘根据本发明的示范性实施例的在入口附近从经重建构数据提取的示范性合并峰值的比较；

图19A描绘根据本发明的示范性实施例的经重建构波形；

图19B描绘根据本发明的示范性实施例的每一频率要素的不加窗位对每一频率要素的加窗位；

图20描绘根据本发明的示范性实施例的渐细窗；

图21描绘根据本发明的示范性实施例的不加窗经重建构波形对加窗经重建构波形；以及

图22为根据本发明的示范性实施例的可经配置为用于执行方法的计算机、系统和/或服务器的计算机的简化功能框图。

具体实施方式

本发明的示范性实施例部分涉及使用化学敏感场效应晶体管（“chemFET”）的大阵列，且更确切地说，涉及离子敏感场效应晶体管（“ISFET”），所述ISFET基于监视在反应期间出现、产生和/或使用的分析物来监视反应，包含例如脱氧核糖核酸（例如，DNA）测序反应。

包含chemFET的大阵列的阵列可用以在多种化学和/或生物过程（例如，生物或化学反应、细胞或组织培养或监视、神经活动、核酸测序等）中检测和测量多种分析物（例如，氢离子、其它离子、非离子分子或化合物等）的静态和/或动态量或浓度，在所述化学和/或生物过程中，可基于此些分析物测量结果获得宝贵的信息。可在检测分析物的方法和/或经由chemFET表面处的电荷改变监视生物或化学过程的方法中使用此些chemFET阵列。因此，本文中论述的系统、方法和计算机可读媒体的至少某些实施例提供chemFET阵列的用途，其涉及检测溶液中的分析物，和/或检测键结到chemFET表面的电荷的改变。

图1描绘根据本发明的示范性实施例的包含大规模chemFET阵列的核酸处理系统。核酸处理系统的实例为核酸测序系统。出于说明的目的，将阵列的chemFET传感器描述为ISFET，所述ISFET是针对对于静态和/或动态离子浓度（包含但不限于氢离子浓度）的敏感度配置的。然而，应了解本发明在这方面不受限制，且在本文中论述的使用ISFET作为说明性实例的实施例中的任一个中，可在替代实施例中类似地使用其它类型的chemFET。类似地，应了解本发明的各种方面和实施例可使用ISFET作为传感器，又检测不是氢离子的一个或多个离子物质。

系统1000可包含半导体/微流体混合结构300，所述混合结构包括ISFET传感器阵列100和微流体流槽200。流槽200可包括安置在ISFET阵列100的对应传感器上方的数个井（未图示）。流槽200可经配置以促进安置于流槽中的一个或多个相同模板核酸的测序，所述测序是经由受控及有序引入到数个测序试剂272（例如，dATP、dCTP、dGTP、dTTP（在本文中一般被称作dNTP）、例如但不限于Mg2+的二价阳离子、冲洗溶液及其类似者的流槽进行的。

如图1中所说明，可经由由计算机260控制的一个或多个阀270和一个或多个泵274来实现将测序试剂引入到流槽200。可使用数种技术将各种处理材料（例如，溶液、样本、反应试剂、冲洗溶液及其类似者）导入（即，引入）到此流槽的井中。如图1中所说明，可将包含dNTP的试剂导入到流槽（例如，经由计算机控制的阀270和泵274），所述试剂从所述流槽扩散到井中，或可通过例如喷墨等其它手段将试剂添加到流槽。在又一实例中，流槽200可不含有任何井，且可利用试剂的扩散性质来限制ISFET阵列100的相应传感器之间的串扰，或核酸可固定在ISFET阵列100的传感器表面上。

图1的系统中的流槽200可以多种方式进行配置以在接近于ISFET阵列100处提供一个或多个分析物（或一个或多个反应溶液）。举例来说，模板核酸可直接附着到或在合适的距离涂覆到传感器阵列100的一个或多个像素，或在位于传感器阵列上方但在反应腔室内的载体材料（例如，一个或多个“珠粒”）中或上，或在自身传感器表面上。处理试剂（例如酶类，例如聚合酶）也可放置在传感器的正上方，或在一个或多个固体载体（例如，所述处理试剂可键结到采集珠粒或其它珠粒）上接近于传感器，或所述处理试剂可成溶液，且自由流动。应理解，可在没有井或珠粒的情况下使用装置。

在图1的系统1000中，根据一个实施例，ISFET传感器阵列100监视离子物质，且确切地说，监视包含氢离子的离子物质的水平/量和/或浓度的改变。物质可由核酸合成或测序反应产生。

本发明的各种实施例可涉及监视/测量技术，所述监视/测量技术涉及ISFET的静态和/或动态响应。应理解，尽管提供核酸合成或测序反应的特定实例以说明chemFET（例如ISFET）的瞬态或动态响应，但如下文所论述，可利用chemFET（例如ISFET）的瞬态或动态响应以监视/感测除了核酸合成或测序反应的特定实例之外的其它类型的化学和/或生物活性。

图2描绘根据本发明的示范性实施例的chemFET传感器阵列的示范性CMOS IC芯片实施方案的框图。如图2中所示，传感器阵列100可耦合到阵列控制器250。可将阵列控制器250制造为“独立”控制器，和/或制造为一个或多个计算机兼容的“卡”，所述卡形成计算机260的部分。阵列控制器250的功能可由计算机260经由接口框252（例如，串行接口、经由USB端口或PCI总线、以太网连接等）控制。

阵列控制器250可向阵列100提供各种供应电压和偏压，以及关于行和列选择、像素输出的取样和数据获取的各种信号。确切地说，阵列控制器250可从阵列100读取包含多路复用相应像素电压信号的一个或多个模拟输出信号（例如，Vout1和Vout2），且接着可将这些相应像素信号数字化以向计算机260提供测量数据，所述计算机又可存储和/或处理数据。在一些实施方案中，阵列控制器250也可经配置以执行或促进各种阵列校准和诊断功能。

如图2中所示，阵列控制器250可向阵列100提供模拟供应电压和接地（VDDA、VSSA）、数字供应电压和接地（VDDD、VSSD）和缓冲输出供应电压和接地（VDDO、VSSO）。在一个示范性实施例中，供应电压VDDA、VDDD和VDDO中的每一个约为3.3伏。在另一实施方案中，供应电压VDDA、VDDD和VDDO可低到约1.8伏。可经由分离传导路径将这些电力供应电压中的每一个提供到阵列100以便于噪声隔离。在另一方面，这些供应电压可源自相应电力供应器/调节器，或这些供应电压中的一个或多个可源自阵列控制器250的电力供应器258中的共同源。电力供应器258也可提供阵列操作所需要的各种偏压（例如，VB1、VB2、VB3、VB4、VBO0、VBODY）和用于阵列诊断和校准的参考电压VREF。

在另一方面，电力供应器258包含一个或多个数/模转换器（DAC），所述DAC可由计算机260控制以允许偏压、参考电压和供应电压中的任一个或全部在软件控制（即，可编程偏置设置）下改变。举例来说，响应于计算机控制（例如，经由软件执行）的电力供应器258可促进调整供应电压中的一个或多个（例如，取决于如由识别码表示的芯片类型在3.3伏与1.8伏之间切换），和/或调整用于像素漏极电流的偏压VB1和VB2、用于列总线驱动的VB3、用于列放大器带宽的VB4和用于列输出缓冲器电流驱动的VBO0中的一个或多个。在一些方面，可调整一个或多个偏压以优化来自启用像素的信号的稳定时间。另外，用于阵列的所有ISFET的共同体电压VBODY可在任选的后制造UV幅射处理期间接地以减少捕获的电荷，且接着在阵列的诊断分析、校准和正常操作期间耦合到较高电压（例如，VDDA）以用于测量/数据获取。同样，可改变参考电压VREF以便于多种诊断和校准功能。

又如图2中所示，通常结合将由阵列100测量的分析物溶液使用的参考电极76（如上文结合图1所论述）可耦合到电力供应器258以提供用于像素输出电压的参考电位。举例来说，在一个实施方案中参考电极76可耦合到供应接地（例如，模拟接地VSSA）以提供用于像素输出电压的参考。在其它示范性实施方案中，参考电极电压可通过以下操作设定：将所关注的具有已知pH水平的溶液/样本放置在接近于传感器阵列100处，以及调整参考电极电压直到阵列输出信号Vout1和Vout2提供处于所要参考电平的像素电压为止，像素电压自所要参考电平起的后续改变反映相对于已知参考pH水平的pH的局部改变。一般来说，应了解与参考电极76相关联的电压不一定要等同于上文所论述的参考电压VREF（其可用于多种阵列诊断和校准功能），但在一些实施方案中，由电力供应器258提供的参考电压VREF可用以设定参考电极76的电压。

关于来自阵列100的数据获取，在一个实施例中，图2的阵列控制器250可包含一个或多个前置放大器253以进一步缓冲来自传感器阵列100的一个或多个输出信号（例如，Vout1和Vout2），且提供可选择的增益。在一方面，阵列控制器250可包含用于每一输出信号的一个前置放大器（例如，用于两个模拟输出信号的两个前置放大器）。在其它方面，前置放大器可经配置以接受从0.0伏到1.8伏或从0.0伏到3.3伏的输入电压，可具有可编程/计算机可选择增益（例如，1、2、5、10和20）和低噪声输出（例如，<10 nV/sqrtHz），且可提供低通滤波（例如，5 MHz和25 MHz的带宽）。关于噪声减少和信噪比的增加，在阵列100经配置为放置在含有阵列控制器250的全部或一部分的印刷电路板的芯片插口中的专用集成电路的一个实施方案中，可在接近于芯片插口处（例如，ZIF插口的底面）使用滤波电容器以便于噪声减少。在又一方面，前置放大器253可具有用于输入和/或输出电压信号的可编程/计算机可选择偏移以将标称电平设定为所希望的范围。

图2的阵列控制器250还包括一个或多个模/数转换器254（ADC）以将传感器阵列输出信号Vout1和Vout2转换成数字输出（例如，10位或12位）以便将数据提供到计算机260。在一方面，一个ADC 254可用于传感器阵列100的每一模拟输出，且每一ADC 254可耦合到对应前置放大器253的输出（在给定实施方案中使用前置放大器的情况下）。在另一方面，ADC254可具有计算机可选择输入范围（例如，50 mV、200 mV、500 mV、1V）以便于与不同范围的阵列输出信号和/或前置放大器参数兼容。在又其它方面，ADC 254的带宽可大于60 MHz，且数据获取/转化率大于25 MHz（例如，高达100 MHz或更大）。

在图2的实施例中，ADC获取定时和阵列行和列选择可由定时发生器256控制。确切地说，定时发生器256提供数字垂直数据和时钟信号（DV、CV）以控制行选择，提供数字水平数据和时钟信号（DH、CH）以控制列选择，且提供列样本和保持信号COL SH以取样启用行的相应像素电压。定时发生器256也将取样时脉信号CS提供到ADC 254以便适当地取样和数字化给定阵列模拟输出信号（例如，Vout1和Vout2）的数据流中的连续像素值。在一些实施方案中，定时发生器256可由执行代码的微处理器实施，且经配置为多通道数字图形发生器以提供适当定时的控制信号。在一个示范性实施方案中，定时发生器256可实施为现场可编程门阵列（“FPGA”）。

如由定时发生器256提供的各种阵列控制信号可用以从传感器阵列100获取像素数据。出于以下论述的目的，“帧”可为包含阵列中的每一像素的值（例如，像素输出信号或电压VS）的数据集，且“帧速率”可为可从阵列获取连续帧的速率。因此，帧速率基本上对应于阵列的每一像素的“像素取样速率”，因为来自任何给定像素的数据是以所述帧速率获得的。

帧速率可为20帧/秒。然而，应了解根据本发明的阵列和阵列控制器在这方面不受限制，因为包含较低帧速率（例如1到10帧/秒）或较高帧速率（例如，25、30、40、50、60、70到100帧/秒等）的不同帧速率是可能的，其中阵列具有相同或较高数目个像素。在一些示范性应用中，可获取包含在若干秒内的许多帧的数据集，以对给定一或多个分析物进行实验。若干此类实验可依次执行，在一些情况下会有暂停，从而允许数据传送/处理和/或清洗传感器阵列以及试剂制备，以供后续实验使用。

举例来说，关于用于检测核苷酸合并的方法，可选择适当的帧速率以充分取样ISFET的输出信号。在一些示范性实施方案中，氢离子信号可具有约1秒到约2.5秒的半峰全宽（FWHM），这取决于核苷酸合并事件的数目。在这些示范性值给定的情况下，20 Hz的帧速率（或像素取样速率）可足以可靠地解析给定像素的输出信号中的信号。再次，主要出于说明的目的提供在此实例中给定的帧速率，且在其它实施方案中可涉及不同帧速率。

关于图2，且如上文所论述，阵列控制器250从阵列100读取包含多路复用相应像素电压信号的一个或多个模拟输出信号（例如，Vout1和Vout2），且接着数字化这些相应的像素信号以将测量数据提供到计算机260。计算机260又可存储和/或处理测量数据。

在实施例中，ADC 254可由定时发生器256经由取样时脉信号CS控制，从而以高数据速率取样输出信号Vout1和Vout2，以提供用于每一像素测量结果的两个或大于两个数字化样本，接着可平均化每一像素测量结果。在实施例中，对于所考虑的每一帧的每一像素，可平均化连续帧中的两个或大于两个像素测量结果。此处，输出为所有考虑的帧的每一像素的平均测量结果。作为此帧平均化技术的结果，可以实现每一像素的噪声的降低。

关于图2，根据本发明的实施例，上文描述的帧平均化技术可出现在阵列控制器250、计算机260或阵列控制器250和计算机260两者中。根据本发明的实施例，计算机260可存储像素测量数据以供进一步处理。在另一实施例中，像素测量数据可存储于在阵列控制器250和计算机260外部的存储器存储装置（未图示）中。

也可对从传感器阵列（例如，图2的阵列100）取样的数据执行可变帧速率平均化。可变帧速率平均化类似于上文所论述的帧平均化技术，还增加了允许可变数目个帧在一起平均化。可变帧速率平均化的益处尤其在于数个帧平均化和输出为单个帧可在数据获取过程期间的任一点进行。一旦对从传感器阵列取样的数据执行可变帧速率平均化，便可对所得数据执行关键帧增量压缩，从而进一步压缩待存储的数据。

随着取样数据速率的加快和传感器阵列密度的升高，像素测量数据可消耗计算机260和/或外部存储器存储装置上的大量存储器。因此期望降低存储器消耗，同时维持像素测量数据的质量。下文详细论述的至少某些示范性方法的目标尤其在于准确地采集与生物/化学事件相关联的数据，同时降低与数据相关联的噪声。此目标可通过实施下文描述的压缩技术来达成。结果，可减少所存储的数据量（例如，在图2的计算机260或外部存储器存储装置中）。

在本发明的一个实施例中，压缩技术可包含处理测序数据且将其存储于频域中。确切地说，频域压缩技术可压缩在空间上相关的井的一小块。举例来说，频域压缩技术可压缩来自约50x50个井的感测阵列的子阵列（“块”）的数据。

可存储截断均值频谱作为关键帧，且可基于关键帧估计每一频率分量的熵。可针对每一频率分量分配数个位以表示个别井与关键帧的差。可计算每一频率分量的同相和正交缩放值，且可按比例缩放井的差。接着可压缩位范围减小的值。

在本发明的一个实施例中，不需要所有频率分量进行适当信号重建构。举例来说，频率分量频谱的前15%到30%可用以产生关键帧，且个别井的频率分量频谱的前15%到30%可用以估计与关键帧的差。

在各种实施例中，个别频率分量可包含比其它频率分量少的唯一信息（例如，具有较低熵）。因此，个别频率分量可需要较少的位来表示其与关键帧的差。在另一实施例中，可舍弃频谱的DC项以移除偏移。在本发明的又一实施例中，可舍弃频域频谱的负频率值，因为源信号表示实数分量。负频率值可被看作正频率分量的复共轭。本发明的至少某些实施例的益处在于本文中论述的压缩技术不需要来自测量之前的时间段的信息或来自中期测量之前的时间段的信息。

在实施例中，每一井的时域波形数据可使用例如傅立叶变换等积分变换转换成其频域表示。图3描绘根据本发明的示范性实施例的时域波形。如图3中所示的原始波形可具有经移除的偏移。波形300描绘与ISFET阵列（例如，图1的阵列100）流体接触的分析物溶液中的一个或多个离子物质的浓度的逐步改变。波形300可表示ISFET阵列对与ISFET阵列流体接触的分析物溶液的离子强度的改变的动态响应。波形300的x轴表示帧号，所述帧号可随时间而变。取决于由定时发生器（例如，图2的定时发生器256）提供的时钟信号，从ISFET阵列对帧进行取样的数据速率可发生变化，如所属领域的一般技术人员应理解。波形300的y轴表示计数数目，计数数目表示由ISFET阵列测量的电压。

图3的波形300可包含ISFET阵列响应“脉冲”，所述脉冲为ISFET阵列特性，所述ISFET阵列特性也被称作“离子步骤或“逐步”响应”。如图3中所示的波形可通过使用积分变换转换成其频域表示。用以将时域波形数据转换成频域波形频谱数据的变换可包含傅立叶变换、傅立叶正弦变换、余弦变换、离散余弦变换、傅立叶余弦变换、哈脱莱变换、梅林变换、双侧拉普拉斯变换、拉普拉斯变换、维尔斯特拉斯变换汉克尔变换、亚伯变换、希尔伯特变换、泊松核和/或恒等变换中的一个或多个。以下文件涉及积分变换，且以全文引用的方式并入本文中：Narasimha, M等人的“关于离散余弦变换的计算（On the Computation of the Discrete Cosine Transform）”，IEEE通信学报，第COM-26卷，第6期，第934页到第936页（1978年）；以及Martucci, S.的“对称卷积和离散正弦和余弦变换（Symmetric Convolution and the Discrete Sine and Cosine Transforms）”，IEEE信号处理学报，第42卷，第5期，第1038页到第1051页（1994年）。

图4描绘根据本发明的示范性实施例的频域波形频谱。如图4中所示的原始波形频谱可通过使用傅立叶变换转换图3的时域波形数据来获得。如图4中所示，实线描绘同相数据，且虚线描绘正交数据。如图4的示范性频域波形频谱可显而易见，频谱的右半边可为频谱的左半边的复共轭（即，正交分量可反相）。

在获得传感器阵列的每一井的频域波形频谱之后，可产生阵列的所有井的平均（“均值”）频谱。图5描绘根据本发明的示范性实施例的传感器阵列中的所有井的平均“均值”频谱。在本发明的一个实施例中，可保留均值频谱的所有频率分量。在本发明的另一实施例中，可保留频率分量的一部分。举例来说，如图5中所示，可保留频率分量的加框部分502。在此实例中，可使用前16的频率分量，且可舍弃其它频率分量（即，17和以上）。在实例实施例中，可使用非DC频率分量，且DC=0 Hz分量可为恒定偏移（图6）。

图6描绘根据本发明的示范性实施例的图5中所示的平均“均值”频谱的截断平均“均值”频谱。截断平均“均值”频谱可为平均频谱的前15%到30%。截断均值频谱可存储为用于压缩技术的关键帧。

一旦已产生关键帧，则可基于关键帧估计个别井的每一频率分量的熵。图7描绘根据本发明的示范性实施例的平均“均值”频谱（关键帧）与个别井之间的频谱相关。接着可识别可生成不良和/或坏数据的井。

在本发明的一个实施例中，来自井的坏数据可通过计算具有复共轭关键帧向量（“K”）的每一井向量（“W”）的相关系数（“C”）来识别，其可由以下公式表示：C = K` * W。如图7中所示，提供良好数据的个别井可形成紧密的聚类。具有低于均值超出一个标准差的“C”值的井可被视为生成坏数据的井。举例来说，在图7的虚线圆702内部的井可低于均值超出一个标准差，且可被舍弃。

具有低于均值超出一个标准差的值的井可作为钉扎、削减和/或异常数据被移除，如图8中所示。图8描绘根据本发明的示范性实施例的基于频谱相关拒绝的井迹线。如图8中所示，顶部线802指示已钉扎的井，中间线804指示提供异常数据的井，且底部线806指示已削减的井。

如上文所提及，对于传感器阵列的每一井，可获得频率要素。接着，对于每一井，频率要素中的每一个可乘以对应关键帧频率要素的复共轭。接着可针对每一相乘频率要素计算群体量值均值、量值标准差和相角标准差。具有更多信息内容（较高熵）的频率要素可具有较大相位和相对量值标准差。

图9A到9I描绘根据本发明的示范性实施例的平均“均值”频谱与个别井之间的各种示范性频谱要素相关。如图9A到9C中所示，其描绘平均“均值”频谱与具有极高熵的个别井之间的示范性频谱要素相关。图9D、9F和9H描绘平均“均值”频谱与具有高熵的个别频谱要素之间的示范性频谱相关。图9G和9I描绘平均“均值”频谱与具有中等熵的个别频谱要素之间的示范性频谱相关。图9E描绘平均“均值”频谱与具有中等到低熵的个别频谱要素之间的示范性频谱相关。

在计算每一频率要素的群体量值均值、量值标准差和相角标准差之后，可计算编码井的频率分量所需要的位数。编码每一频率要素增量值（最小所需位[minBits]到最大所需位[maxBits]）所需要的位数可基于每一频率要素的熵值来计算。每一频率要素n的熵值（“

”）可等于频率要素n的相角标准差（“

”）乘以频率要素n的标准差的量值（“

”），如以下公式所描绘：

。

在计算每一频率要素的熵值之后，可归一化熵值（即，经缩放以具有最大值1）。接着，可将归一化熵值转换成位。每一频率要素n所需要的位（“

”）可等于（所需要的最大位减去1（“

”)）加频率要素n的归一化熵值的log 2（“

”），如由以下公式所描绘：

。bits _n的最小值可限于（minBits-1）。另外，考虑到位值的正负号，每一频率要素n所需要的位可增加一（1）。

图10描绘根据本发明的示范性实施例的归一化频率分量熵。图11描绘根据本发明的示范性实施例的每一频率分量的位。在实例实施例中，minBits可等于三（3），且maxBits可等于九（9）。此实例的数据可经编码为带正负号的整数，这是计算和使用（minBits-1）和（maxBits-1）限制量值的原因，其中在结尾处添加正负号位。当然，所属领域的一般技术人员应理解可向minBits和maxBits指派其它值。

在计算关键帧的情况下，可计算每一井的截断频谱与关键帧频谱之间的差。图12A和12B分别描绘根据本发明的示范性实施例的示范性关键帧和井的截断频谱。图13描绘根据本发明的示范性实施例的来自关键帧的示范性井频谱增量。

在计算与关键帧的井频谱差之后，可计算同相数据的缩放向量（“

”），以及正交数据（“

”）频率分量增量值的缩放向量。首先可确定每一频率分量的同相数据的所有增量的最大幅值（“

”）。接着，可计算每一频率分量n的同相数据的尺度。每一频率分量的同相数据的所有增量的最大幅值（“

”）可除以每一频率要素n所需要的位（“

”）减一的二次幂再减一，其可由以下公式表示：

。

首先可确定每一频率分量的正交数据的所有增量的最大幅值（“

”）。接着，可计算每一频率分量n的正交数据的尺度。每一频率分量的正交数据的所有增量的最大幅值（“

”）可除以每一频率要素n所需要的位（“

”）减一的二次幂再减一，其可由以下公式表示：

。

在上文论述的实例中，同相数据和正交数据的缩放向量可使用bits _n-1来计算，因为计算的值可为无符号量值。另外，bits _n-1可用以确保带正负号的值预留了一位。

在计算同相数据和正交数据频率分量增量值的缩放向量之后，可编码频率增量值。为了编码频率分量的同相数据（“

”），每一频率分量同相数据值（“

”）可除以同相数据的缩放向量（“

”），且接着可舍入到最近整数，如以下公式中所示：

。为了编码频率分量的正交数据（“

”），每一频率分量正交数据值（“

”）可除以正交数据的缩放向量（“

”），且接着可舍入到最近整数，如以下公式中所示：

。

接着可封装经编码频率增量值。来自每一经编码值的推荐数目个经编码位可封装到常规数据字中以供存储。可通过使用位移位和/或通过使用逻辑或来执行封装。举例来说，封装可如下表示：Data2write = encodedI[0] | (encodedQ[0] << bits_n[0]) |encodedI[1] | (encodedQ[1] << bits_n[1]) |等等。

封装数据可存储有标头，且传感器阵列的每一子阵列可具有相应标头和封装数据。标头可包含每一井的频率要素的数目、原始时域样本的数目、关键帧数据（对于同相数据和正交数据两者）、缩放向量（对于同相数据和正交数据两者）以及位/要素向量。在一个实施例中，可顺序封装和存储经压缩的井数据。可顺序存储给定井的位，且井可按顺序光栅（或任何预定义）次序存储。每一经压缩井所需要的位数可将“位/要素”向量之和加倍。

图14描绘根据本发明的示范性实施例的压缩方法1400。在步骤1402中，可接收与每一井或传感器阵列的其它反应区域中发生的化学事件相关联的模拟输出或波形。在步骤1404中，可使用数字转换器将模拟输出或波形转换成包括多个帧的数字化输出或波形。

数字转换器可具有任何合适的数字转换器类型，且可例如将模拟输出或波形转换成14位数据。可使用各种转换帧速率。举例来说，速率可为0.008478秒/帧（例如，IonTorrent™的离子314™芯片）、0.034402秒/帧（例如，Ion Torrent™的离子316™芯片）、0.061075秒/帧（例如，Ion Torrent™的离子318™芯片）以及0.033秒/帧（例如，IonTorrent™的离子质子™芯片I或离子P1™芯片。更一般地说，速率可在约0.004秒/帧与约0.333秒/帧之间选择。

在替代实施例中，可接收包括多个帧的数字化输出或波形数据。波形数据可包含具有多个帧的每一井的基于时域的波形。

在步骤1406中，可将每一井的数字化的基于时域的波形转换成频域谱。可通过使用例如傅立叶变换等积分变换将每一井的数字化的基于时域的波形转换成频域谱。接着，在步骤1408，关键帧可从来自所有井的频域谱数据产生。在一个实施例中，关键帧可为所有井的平均（均值）频谱。或者，平均“均值”频谱的截断平均“均值”频谱可用以产生关键帧。举例来说，可从平均频谱的前15%到30%产生关键帧。

在步骤1410中，对于传感器阵列的每一井，可计算关键帧与井之间的差。在计算关键帧与井中的每一个之间的差之后，在步骤1412，可编码所述差，如上文所论述。

下文提供用于压缩测序数据的源代码的实例：

function [keyframe, bitsPerFreq, scale, deltas, badIdx] = DFCompress(rawData, freqPts, maxBits)

%步骤1-转换成频域

tmp = zeros (size(rawData)); %预分配

for idx = 1:size(rawData,2)

tmp(:,idx) = fft(rawData(:,idx));

end

%步骤2-舍弃DC和负频率分量，因为DC和

%负频率分量仅为正%频率的复共轭。仅保留前指定数目个频率样本。

freqData = tmp(2:(freqPts+1),:);

%步骤3-产生关键帧和相关统计

keyframe = zeros(1,size(freqData,1)); %预分配

for idx = 1 :3izc(frcqDntn, l)

keyframe(idx) = mean(freqData(idx,:));

end

corrKnl = keyframe'; %使用复共轭作为相关核

corrVal = zeros(1,size(freqData,2)); %预分配

z = zeros(size(freqData)); %预分配

for idx = 1:size(freqData,2)

corrVal(idx) = freqData(:,idx).' * corrKnl; %相关井

z(:,idx) = freqData(:,idx) .* corrKnl; %相关频率分量

end

corrNorm = corrVal ./ (keyframe*corrKnl); %归一化

badidx = find(abs(corrNorm) < (1-std(abs(corrNorm)))); %识别坏井

avgMag = zeros(1,size(freqData,1)); %预分配

sdMag = zeros(1,size(freqData,1)); %预分配

sdAng = zeros(1,size(freqData,1)); %预分配

for idx = 1:size(freqData,1)

avgMag(idx) = mean(abs(z(idx,:)));

sdMag(idx) = std(abs(z(idx,:)));

sdAng(idx) = std(angle(z(idx,:)));

end

%步骤4-估计每一频率点的熵且转换成位/

%样本（数据为复数，因为总位数将是双倍）

%总位数为正负号+(maxBits-1) 2的互补量值

%前三个样本通常将为零，因为（一般来说）它们

%仅含有关于背景pH级的信息。

rawEmphasis = (sdAng .* sdMag) ./ avgMag;

emphasis= rawEmphasis ./ max(rawEmphasis); %归一化

bitsPerFreq = round((maxBits-1)+log2(emphasis)); %每一值的量值位

bitsPerFreq(find(bitsPerFreq < 2)) = 2 ; %约束最小位数

bitsPerFreq = bitsPerFreq + 1 ; %添加正负号位

%步骤5-计算非零bitsPerFreq频率的原始频谱增量值和

%缩放参数

encldx = 1;

bitsPerFreqIdx = find(bitsPerFreq) ;

rawDelta = zeros(length(bitsPerFreqIdx),size(freqData,2)); %预分配

deltaRange = zeros(l,length(bitsPerFreqIdx)); %预分配

scale= zeros(1,length(bitsPerFreqIdx)); %预分配

for idx = bitsPerFreqIdx

rawDelta(encIdx,:) = freqData(idx,:) - keyframe(idx);

tmpReal = max(abs([max(real(rawDelta(encIdx,:))),

min(real(rawDelta(encIdx,:)))]));

tmpImag = max(abs([max(imag(rawDelta(encldx,:))),

min(imng(rawDelta(encldx,:)))]));

deltaRange(encIdx) = tmpReal + 1i*tmplmag;

scale(encIdx) = deltaRange(encIdx)/((2^(bitsPerFreq(idx)-1))-1);

encldx = encIdx + 1;

end

%步骤6-量化处于最小/最大界限处的频谱增量值和削减值

deltas = zeros(size(rawDelta)); %预分配

for idx = 1:size(rawDelta,1)

deltas (idx,:) =

round((real(rawDelta(idx,:))/real(scale(idx)))+1i*(imag(rawDelta(idx,:))/

imag(scale(idx)))) ;

end

%步骤7-基于bitsPerFreq进行位封装。

在已编码和存储井数据之后，关键帧和增量信息可用以重建构每一井的频谱。为了准确重建构每一井的频谱，可从正频谱的复共轭重建构负频率。此外，高阶频率分量（正和负）可用零填充。在已重建构频域谱之后，积分变换的逆可用以将经重建构频谱转换回到时域波形。举例来说，傅立叶逆变换可用以将经重建构频谱转换回到时域波形。经重建构时域波形可为纯实数零均值（无复分量）。

下文提供用于解压缩（扩展）测序数据的源代码的实例：

%差分傅立叶压缩-扩展器

function [data] = DFExpand(keyframe, bitsPerFreq, scale, deltas,timePts)

%步骤1-建构频谱模板

freqPts = length(bitsPerFreq);

posLength = ceil(timePts / 2);

templateFreq = zeros(1, posLength);

templateFreq(2:(freqPts+1)) = keyframe;

%步骤2-将增量值添加到模板以重建构井频谱的

%正半边

numWells = size(deltas,2);

scaleIdx = find(bitsPerFreq)+1; %得到增量的频率区间索引值

freqData = zeros(posLength, numWells); %预分配

for idx = 1:numWells

freqData(:,idx) = templateFreq;

freqData(scaleIdx,idx) = freqData(scaleIdx,idx) + (real(scale.’) .*

real(deltas(: idx))) + 1i*(imag(scale.’) * imag(deltas(:,idx)));

end

%步骤3-从正频谱数据的复共轭合成

%负频谱数据

freqData((posLength+1):timePts,:) = (freqData(posLength:-1:2,:)’).’;

data= zeros(size(freqData)) ; %预分配

for idx = 1:numWells

data(:,idx) = ifft(freqData(:,idx));

end

图15描绘根据本发明的示范性实施例的经重建构井频谱。对于传感器阵列中的每一井，可重建构井频谱。可将经重建构频谱样本的总数设定为等于时间样本的数目。接着，可根据增量值乘以缩放值来建构正低频复频谱，其中添加了关键帧。DC和正高频值可为零，且负频率分量可从正频谱的镜像复共轭建构。

图16描绘根据本发明的示范性实施例的经重建构时域波形。在建构传感器阵列中的井中的每一个的井频谱之后，可使用逆积分变换将频域谱转换成时域波形。举例来说，傅立叶逆变换可将经重建构频域谱转换回到时域波形。如上文所提及，时域波形可为纯实数，且具有零均值。可通过在应用傅立叶逆变换之前将DC频率要素设定为非零值而引入偏移。DC频率值可为偏移除以样本数目，且可为纯实数（无复分量）。

在实例实施例中，来自离子质子™系统的开发用于下一代测序的技术原型芯片的处于105个样本/流的ISFET的数据在压缩到25个频率要素之后经重建构，其中跨越440个流，minBits等于三（3），maxBits等于九（9）。传感器阵列的井的50乘50群组实现每一流的每一井的27.69字节的平均压缩。图17A到17J描绘根据本发明的实施例的传感器阵列的各种井的示范性原始时域波和经重建构时域波形。如图17A到17J中可以看出，在大多数情况下，原始波形和经重建构波形几乎不可区分。

下文提供用于测试测序数据的压缩的压缩比的源代码的实例：

%用以测试单个流的DFC压缩算法的函数

function [data, raw, compressionRatio, bytesPerWell] =

testFlowDFC (fileName, patchWidth , patchHeight, maxBits, freqPts ,wellUdx)

idxVal = 11:60;

%上传缩略图数据

[img,ftimes]=Loadimage(fileName, [0 0 patchHeight patchWidth]);

analysisRegion=img(1:patchHeight,1:patchWidth,:); % extract 50x50 block

%压缩

timePts = size(anal ysisRegion, 3); %预期每一井105个帧

rawData=reshape(analysisRegion,patchWidth*patchHeight,timePts).’;

%压缩数据

[keyframe, bitsPerFreq, scale, deltas, badIdx] = DFCompress(rawData,

freqPts, maxBits);

%计算压缩比

rawDataSize = size(rawData,1) * size(rawData,2) * 2; % unpacked bytes

keyframeSize = freqPts * 4 ; % 16位I+16位Q

bitsPerFreqSize = freqPts * 0.5; %每一频率点4位半字节

scaleVecSize = length(scale) * 4; %每一非零bitsPerFreq点的16位I+16位Q

%缩放

deltaSize = (sum(bitsPerFreq) / 4) * size(rawData,2); %封装字节

compDataSize = keyframeSize + bitsPerFreqSize + scaleVecSize +deltaSize;

compressionRatio = rawDataSize / compDataSize ;

bytesPerWell = compDataSize / size(rawData,2);

%重建构数据

[reconData] = DFExpand(keyframe , bitsPerFreq, scale, deltas ,timePts);

%产生移除偏移的原始数据以用于比较

offsetData = zeros(size(rawData)); %预分配

for idx = 1:size(rawData,2)

tmp = fft(rawData(:,idx));

tmp(1) = 0 ; %移除DC偏移

offsetData(:,idx) = ifft(tmp);

end

%提取数据以供进一步分析

data = reconData(idxVal, wellIdx);

raw = offsetData(idxVal , wellIdx);

下文提供用于测试测序数据的压缩的压缩比的源代码的另一实例：

%满载测试DFC压缩算法，提取所有所要

%统计。

function [data, raw, compressionRatio, bytesPerWell] =

testFullDFC(dataPath, flows , patchWidth, patchHeight, maxBits,freqPts ,

wellIdx, showMovie)

%预分配缓冲器

data = zeros(50, flows);

raw = zeros(50, flows);

compressionRatio = zeros(1, flows);

bytesPerWell = zeros(1, flows);

h = waitbar(0,['Processing flow 0/’,num2str(flows)]);

for idx = 1:flows

fileName = [dataPath, sprintf(‘/acq_%04d.dat’,(idx-1))];

[data(:,idx), raw(:,idx), compressionRatio(idx),bytesPerWell(idx)] =

testFlowDFC(fileName, patchWidth, patchHeight, maxBits, freqPts,wellldx);

waitbar(idx/flows, h, [‘Processing flow’, num2str(idx),’/’,num2str(flows)]);

end

close(h);

fprintf('%.2f bytes per well\n’, mean(bytesPerWell));

% plot traces

ymax=max([max(max(abs(raw)));max(max(abs(data)))])*1.1;

h = figure(‘position’,[0 0 2000 1600]);

set(h,’PaperUnits’,’inches’,’PaperPosition’, [0 0 20 16])

for figIdx=1:(flows/20)

for y=1:4

for x=1:5

subplot(4,5,x+(y-1)*5) ;

plot(11:60, raw(:,x+(y-1)*5+(figIdx-1)*20), ‘r’); hold on;

plot(11:60, data(:,x+(y-1)*5+(figIdx-1)*20), ‘b’); hold off;

legend(‘Raw’,’DFC’,’location’,’southeast’);

xlim([10 60]); ylim((-ymax ymax]);

title([‘Flow ‘, num2str(x+(y-1)*5+(figIdx-1)*20)]);

end

drawnow;

fname = sprintf(‘Well_%d_%dpts_%dbit_flows_%d_to_%d.png’,wellIdx,freqPts,maxBits,

(figIdx-1)*20+1,figIdx*20);

print(h, ‘-dpng’, ‘-r100’, fname);

end

close (h) ;

% generate “movie”

if exist(‘showMovie’, ‘var’)

h = figure();

for idx=1:flows

figure(h);

plot(11:60,[data(:,idx),raw(:,idx)]);

ylim([-ymax ymax]);

title([‘Well ‘, num2str(wellIdx), ‘, ‘, num2str(freqPts), ‘pts, ‘,

num2str(maxBits), ‘bit limit (flow ‘, num2str(idx), ‘)’]);

drawnow;

pause(0.5);

end

图18A和18B描绘根据本发明的示范性实施例的在入口附近从经重建构数据提取的示范性合并峰值的比较。如图18A和18B中所示，原始波形与用不同频率要素压缩的经重建构波形之间的差可能较大，其中上升时间最尖锐，且在较高频率处具有相当多的信息内容。也如图18A和18B中所示，更精确重建构的波形在与原始波形相比较时可在使用更多频率分量时生成，其又减小压缩比。

根据本发明的实施例，高质量测序可需要约20%到25%的频率要素（例如，与未经压缩数据相比较，小于1%的吞吐量下降）。另外，在本发明的实施例中，maxBits值可为至少7位，这可确保用于准确重建构的足够动态范围。在本发明的一个实施例中，编码105个帧的数据可通过具有每一样本最大8位的25个频率要素来实现。因此，可以实现平均跨越440个流的（约来自16位原始的7.5x）压缩每一流的每一井的221.5个位，其中AQ17或AQ20测序性能损失最小。

图19A描绘根据本发明的示范性实施例的经重建构波形。如图19A中所示，经重建构波形可具有“摆动”1902，其是由原始波形不在零值处开始引起。摆动可将系统残差误差引入到经重建构波形数据中。图19B描绘根据本发明的示范性实施例的每一频率要素的不加窗位对每一频率要素的加窗位。在本发明的一个实施例中，在使用例如傅立叶变换等积分变换之前，渐细窗可作为缩放向量应用到DC移除原始波形。渐细窗的逆接着可在例如傅立叶逆变换等逆积分变换之后应用于经重建构数据。

图20描绘根据本发明的示范性实施例的渐细窗。如图20中所示，渐细窗可包含三个部分：前斜坡2002、中间部分2004和尾斜坡2006。在实例实施例中，前斜坡2002和尾斜坡2006可为均匀长度高斯窗的一半。斜坡的长度和曲率可以是可配置的。在一个实施例中，曲率可由高斯参数“阿尔法”控制。在实例实施例中，中间部分可为“整体增益”，其中可不对原始波形数据作出改变。

图21描绘根据本发明的示范性实施例的不加窗经重建构波形对加窗经重建构波形。如图21中所示，可经由使用渐细窗显著地减少摆动。在与不加窗数据相比较时，可改进压缩以获得等效测序性能。较短窗长度也可实现较佳结果。在一个实例实施例中，前五个（5）帧和上五个（5）帧用于斜变过渡。在实例实施例中，最优“阿尔法”参数可为约2.0+/-15%到约2.3+/-15%（用以产生前/尾斜坡的高斯窗等式）。在特定实施例中，最优“阿尔法”参数可为约2.15+/-15%。在示范性实施例中，AQ47性能相比不加窗结果可得到改进。

下文提供用于使用测序数据的压缩的渐细窗的源代码的实例：

#include <iostream>

#include <algorithm>

#include “DfcCompr.h”

#include <malloc.h>

#include <cmath>

using namespace std;

/***

DFT引擎可保持与帧数的任何改变同步。

*/

int DfcCompr::SetNumFrames(int frames)

/***

if ((frames> 0) && (frames != n_frames))

{

n_ frames = frames;

DFT.SetPts(frames);

}

return n_frames;

}

/***

范围检查新窗长度值且产生新窗

系数。窗长度可能相对较短。

*/

int DfcCompr::SetWindowLength(int length)

{

if ((length >= 0) && (length < (n_frames/4)) && (length !=winParamLen))

{

winParamLen = length;

GenerateWindow();

}

Return winParamLen;

}

/***

范围检查新窗阿尔法值且产生新窗系数。

*/

float DfcCompr::SetWindowAlpha(float alpha)

{

if ((alpha >= 0.0f) && (alpha <= 3.0f) && (alpha !=winParamAlpha))

{

winParamAlpha = alpha;

GenerateWindow();

}

return winParamAlpha;

}

/***

基于当前窗窗度集合和阿尔法配置参数产生

高斯窗系数的新集合。

*/

void DfcCompr::GenerateWindow()

{

//设置向量以存储压缩和解压缩窗

windowCompress.resize(winParamLen*2);

windowExpand.resize(winParamLen*2);

float halfN = static_cast<float>(winParamLen);

float n = 0.5f-halfN; //在-(N-1)/2开始

int idx, maxidx;

float an_ sq, val, oneOnVal;

maxidx = (winParamLen * 2) - 1 ;

for (idx = 0; idx < winParamLen; ++idx)

{

an_sq = winParamAlpha * (n / halfN);

val = exp((-0.5f) * (an_sq * an_sq));

oneOnVal = 1.0f / val;

n += 1.0f; // 为下一循环迭代设置n

windowCompress[idx] = val;

windowCompress[maxIdx- idx) = val; //利用窗对称性

windowExpand [idx] = oneOnVal ;

windowExpand [maxIdx-idx] = oneOnVol;

}

/***

将图像数据分解成基向量的集合且每一井投影到其上。

@param n_wells –图像分块中的井数目

@par am n_frame -图像分块中的帧数。

@par am image –与Rawlmage结构一样的次序帧中的数据的个别井，row、col主次序，因此value(row_ i,col_ j,frame_k) = image[row_ i * ncol + col_j+ (nrow * ncol * frame_k)]

@param n_sample_wells –分块的样本中的井数目

@param image_sample –向量上方的图像的样本，且如此

@param compressed –输出有损失的经压缩分块

*/

void DfcCompr::LossyCompress(float *image)

{

//////////////////////////////////////////////////

//临时scratch_pad存储器分配

if (scratch_pad == NULL) scratch_pad = new float[n_wells*n_basis*2];

//////////////////////////////////////////////////

int wellldx, frameldx;

float tmpWell[n_frames];

register float tmpSum;

register int offset, winldx;

float * __restrict dftPtrI = &scratch_pad[0];

float * __restrict dftPtrQ = &scratch_pad[n*wells*n_basis];

float deltaMagI[n_basis);

float deltaMagQ(n_basis);

// Initialize keyframe

for (frameIdx = 0; frameIdx < n_basis; ++frameIdx)

{

keyFrameI[frameIdx] = 0.0f;

keyFrameQ[frameIdx] = 0.0f;

}

for (wellIdx = 0 , offset = 0 ; wellIdx < n_wells; ++wellIdx)

{

//步骤1a–从图像立方体提取下一井

tmpSum = 0.0f;

winIdx = 0;

for (frameIdx = 0 ; frameIdx < n_frames; ++frameIdx)

{

tmpWell[frameIdx] = image[(frameIdx*n_wells)+wellIdx];

tmpSum += image[(frameIdx*n_wells)+wellIdx);

}

//步骤1b-移除DC偏移且将加窗应用于每一井

tmpSum /= static_cast<float>(n_frames); //DC偏移

//应用窗的前部

for (frameIdx = 0 ; frameIdx < winParamLen ; ++frameIdx)

{

tmpWell[frameIdx] -= tmpSum;

tmpWell[frameIdx) *= windowCompress[winIdx++];

}

//应用中间部分（窗是一整体)

for (; frameIdx < (n_frames-winParamLen); ++frameIdx)

tmpWell[frameIdx] -= tmpSum;

//应用窗的后部

for (; frameIdx < n_ frames; ++frameIdx)

{

tmpWell[frameIdx) -= tmpSum;

tmpWell[frameIdx) *= windowCompress[winIdx++];

//步骤2-计算部分DFT

//存储增量且累计关键帧

winIdx = n_basis * wellIdx;

DFT.PartialDFT(1, static_ cast<unsigned int>(n_basis),

&tmpWell[0], &dftPtrI[winIdx], &dftPtrQ[winIdx]);

for (frameIdx = 0; frameldx < n_basis; ++frameIdx, ++offset)

{

keyFrameI[frameIdx] += dftPtrI[offset];

keyFrameQ[frameIdx] += dftPtrQ[offset];

}

//步骤3-将累计的频谱转换成均值关键帧，且

//初始化

//最大增量量值向量

for (frameIdx = 0 ; frameIdx < n_basis; ++frameIdx)

{

keyFrameI[frameIdx] /= static_cast<float>(n_wells);

keyFrameQ[frameIdx] /= static_cast<float>(n_wells);

deltaMagI[frameIdx] = 0.0f;

deltaMagQ[frameIdx] = 0.0f;

}

//步骤4-计算相关统计和所需要的重要向量

//以充填每一频率要素向量的位。

Emphasis();

//步骤5-将DFC数据从绝对值转换成原始增量，一直保持

//每一频率分量所观测到的

//最大分量量值

for (wellIdx = 0, offset = 0; wellIdx < n_ wells; ++wellIdx}

//能够在此处使用SSE/AVX向量加速来并行处理

//多个频率要素

for (frameIdx = 0; frameIdx < n_ basis; ++frameIdx, ++offset)

{

dftPtrI[offset] -= keyFrameI[frameIdx];

dftPtrQ[offset] -= keyFrameQ[frameIdx];

if (deltaMagI[frameIdx] < abs(dftPtrI[offset]))

deltaMagI[frameIdx] = abs(dftPtrI[offset]);

//上传最大观测到的实数量值向量

if (deltaMagQ[frameIdx ] < abs(dftPtrQ[offset]))

deltaMagQ[frameIdx] = abs(dftPtrQ[offset]);

//上传最大观测到的图像量值

//向量

}

//充填缩放向量

for (frameIdx = 0 ; frameIdx < n_basis; ++frameIdx)

{

tmpSum = static_cast<float>((1<<static_cast<int>(bitsPerFreq[frameIdx]-1))-1);

scaleVectorI[frameIdx] = deltaMagI[frameIdx] /tmpSum;

scaleVectorQ[frameIdx] = deltaMagQ[frameIdx] /tmpSum;

//步骤6-量化频谱增量值且封装至输出

//向量中

for (wellIdx = 0, offset = 0; wellIdx < n_ wells; ++wellIdx)

//能够在此处使用SSE/ AVX向量加速来并行处理

//多个频率要素

for (frameIdx = 0; frameIdx < n_basis; ++frameIdx, ++offset)

{

deltaI[offset] = static_cast<short>(round(dftPtrI[offset]/scaleVectorI[frameIdx]));

deltaQ[offset] = static_cast<short>(round(dftPtrQ[offset]/scaleVectorQ[frameIdx]));

}

//////////////////////////////////////////////////

//临时scratch_pad存储器解分配

delete scratch_pad; scratch_pad = NULL;

//////////////////////////////////////////////////

}

void DfcCompr::Emphasis()

{

float kernel[n_basis];

float *dftPtrI = &scratch_pad[0];

float *dftPtrQ = &scratch_pad[n_wells*n_basis];

float *corrI =new float[n_wells*n_basis];

//可使其为静态或半静态

float *corrQ = new float[n_wells*n_basis];

// 可使其为静态或半静态

register int offset, freqIdx;

int wellIdx;

float meanMag[n_basis];

float meanAng[n_basis];

float varMag[n_basis];

float varAng[n_basis];

register float delta;

float emphasisVector [n_basis];

register int bitVal;

//从关键帧产生相关核，且初始化统计

//向量

for (freqIdx = 0; freqIdx < n_basis; ++freqIdx)

{

kernelI[freqIdx] = keyFrameI[freqIdx];

kernelQ[freqIdx] = - keyFrameQ[freqIdx];

meanMag[freqIdx] = 0.0f;

meanAng[freqIdx] = 0.0f;

varMag[freqIdx] = 0.0f;

varAng[freqIdx] = 0.0f;

}

//使每一井频谱与均值井频谱的复共轭

//相关

for (wellIdx = 0, offset = 0 ; wellIdx < n_ wells; ++wellIdx)

//可能够向量化相关过程

for (freqIdx = 0 ; freqIdx < n_basis; ++freqIdx,++offset)

{

// DFT data ==> (a+ib)

// kernel ==> (x+ib)

// correlated data ==> (a+ib) (x+iy) =ax+iay+ibx-by = (ax-by)+i(ay+bx)

corrI[offset] = dftPtrI[offset]*kernelI[freqIdx]-dftPtrQ[offset]*kernel[freqIdx]; // (ax- by)

corrQ[offset] = dftPtrI[offset]*kernelQ[freqIdx]+dftPtrQ[offset]*kernelI[freqIdx]; // (ay+bx)

}

//找到每一频率要素的量值和角度的均值

for (wellIdx = 0, offset = 0 ; wellIdx < n_wells; ++wellIdx)

//可能够向量化频率要素

for (freqldx = 0 ; freqldx < n_basis; ++freqldx,++offset)

{

meanMag[freqIdx] += sqrt((corrI[offset]*corrI[offset])+(corrQ[offset]*corrQ[offset]));

meanAng[freqIdx] += atan2(corrQ[offset],corrI[offset]);

}

for (freqldx = 0; freqldx < n_basis; ++freqldx)

{

meanMag[freqIdx] /= static_cast<float>(n_wells);

meanAng[freqIdx] /= static_cast<float>(n_wells);

}

//找到标准差

for (wellldx = 0, offset = 0; wellIdx < n_ wells; ++wellldx)

//可能够向量化频率要素

for (freqldx = 0; freqIdx < n_basis; ++freqIdx, ++offset)

{

delta = sqrt((corrI[offse]*corrI[offse])+(corrQ[offset]*corrQ[offset]))-

meanMag[freqIdx];

varMag[freqIdx] += delta * delta;

//累计量值方差

delta = atan2(corrQ[offset], corrI[offset])-meanAng[freqIdx];

varAng[freqIdx] += delta * delta;

//累计角度方差

}

//结束相关数据（一旦使此存储器块为静态/半静态，

//便不需要释放）

delete corrI;

delete corrQ;

//将累计的方差转换成原始重要（原始Matlab

//计算使用标准差而非方差，但

//使用此处已经优化的额外sqrt步骤）

for (freqIdx = 0, delta = 0.0f; freqIdx < n_basis; ++freqIdx)

{

emphasisVector[freqIdx] = sqrt((varMag[freqIdx]*varAng[freqIdx])/static_cast<float>(n_wells));

if(emphasisVector[freqIdx] > delta)

delta= emphasisVector[freqIdx];

//保留归一化的最大值

}

//归一化重要向量，且将其转换成每一频率要素的

//位

for (freqIdx = 0; freqIdx < n_basis; ++freqldx)

{

//计算量值位数以进行分配（log将

//总是传回零或负值）

bitVal = static_cast<int>(static_cast<float>(n_maxBits)+ceil(log2(emphasisVector[freqldx]/delta)) ;

if (bitVal < n_minBits) bitVal = n_minBits;

//确保分配至少最小数目个量值

//位

++bitVal ; //添加正负号位

bitsPerFreq[freqIdx] = static_cast<unsigned char>(bitVal);

}

/***

从经压缩的频域向量重建构图像

数据。

@param n_wells –图像分块中的井数目。

@param n frame –图像分块中的帧数。

@param image –与Rawlmage结构相同的次序帧中的数据

的个别井，row、col主次序，因此value(row_i,col_j,frame_k) =

image[row_ I * ncol + col_j + (nrow * ncol * frame_k)]

*/

void DfcCompr::LossyUncompress(float *image)

{

float wellBufI[n_basis];

float wellBufQ [n_basis];

float tmpBuf[n_frames];

int wellldx, frameIdx;

register int offset, winldx;

for (wellldx = 0, offset = 0 ; wellldx < n_wells; ++wellldx)

{

//从关键帧、增量和缩放信息重建构部分

//频谱

for (frameldx = 0; frameldx < n_basis; ++frameldx, ++offset)

{

wellBufI[frameIdx) = keyFrameI[frameIdx]+ scaleVectorI[frameIdx]*static_cast<float>(deltaI[offset]);

wellBufQ[frameIdx] = keyFrameQ[frameldx)+ scaleVectorQ[frameldx]*static_cast<float>(deltaQ[offset]);

}

//使用IDFT从部分频谱重建构时域信号

DFT.PartialIDFT(1, n_basis, &tmpBuf [0], &wellBufI[0], &wellBufQ[[0]);

//将经重建构的值存储于图像立方体中

for (frameIdx = 0, winldx = 0; frameldx <winParamLen; ++frameldx)

image[(frameldx*n_wells)+wellIdx) =tmpBuf[frameIdx] * windowExpand[winldx++];

for (; frameldx < (n_frames-winParamLen); ++frameldx)

image[(frameIdx*n_wells)+wellIdx) =tmpBuf[frameIdx);

for (;frameldx < n_frames; ++frameldx)

image[(frameIdx*n_wells)+wellIdx) =tmpBuf[frameIdx) * windowExpand[winIdx++];

}

#ifndef DFCCOMPR_H

#define DFCCOMPR_H

#include <stdio.h>

#include <string>

#include <vector>

#include "ParallelDFT.h"

/***

由压缩产生的数据可用以解压缩所得数据

*/

class DfcCompr {

public:

DfcCompr() :DFT() {

n_wells = n_ rames = n_ basis = 0; n_maxBits = 7;n_minBits = 2;

basis_vectors = coefficients = scratch_pad =NULL;

keyFrameI = keyFrameQ = scaleVectorI =scaleVectorQ = NULL;

bitsPerFreq = NULL;

deltaI = deltaQ = NULL;

winParamLen = 5; winParamAlpha = 2.15f;

GenerateWindow();

}

~DfcCompr()

{

windowCompress.clear();

windowExpand.clear();

}

/***

n_frames必须隐藏，因为

DFC引擎内部的其它事物现在依赖于n_frames，且还需要在其改变值时

进行重新配置。

*/

int SetNumFrames(int frames);

int n_wells; /// <区中的列数

int n_basis ; /// <偏置向量的数目（频率要素），必须小于帧数的一半

int n_maxBits; /// <每一频率要素的最大位数（排除正负号位）

int n_minBits; /// <每一频率要素的最小位数（排除正负号位）

/***

按串列形式的基向量的矩阵。第一nframe

浮点是第一基向量，接下来的nframe浮点是第二

基向量，等等。

*/

float *basis_vectors; // NOT USED BY DFC

/***

每一基向量的系数矩阵。概念上，

列主格式的nbasis向量的井数目的矩阵（nrow * ncol）

（第一nrow *ncol是第一基向量的

系数）。出于历史原因，井自身

呈行主格式，因此井在位置row_i、col_j处的第n个基向量的

系数（basis_n）将为：

value(row_ i,col_j,basis_n) = \

coefficents[row_i * ncol + col_ j + (nrow * ncol * basis_n)] ;

*/

float *coefficients; // NOT USED BY DFC

/***

针对压缩的井的每一块存储DFC项

- keyFrameI : float * n_basis，vect或含有均值频谱的实数分量

分量

- keyFrameQ : float * n_basis，向量含有均值频谱的虚数

分量

- bitsPerFreq : uchar * n_basis，向量含有用以编码每一频率要素增量

的位数

- scaleVectorI : float * n_basis，向量含有应用于实数频率增量的

共同缩放值

- scaleVectorQ : float * n_basis，向量含有应用于虚数频率增量的

共同缩放值

*/

float *keyFrameI;

float *keyFrameQ;

unsigned char *bitsPerFreq;

float *scaleVectorI;

float *scaleVectorQ;

/***

用于位封装和存储的DFC主体数据。组织每一向量

为一个井的n_basis样本，后跟着

下一井的n_basis样本，等等。

deltaI[n]为“实数”部分，且deltaQ[n]为复增量样本“delta[n]”的

“虚数”部分。

- deltaI : short * n_ basis * n_wells

- deltaQ : short * n_ basis * n_wells

*/

short *deltaI;

short *deltaQ;

/***

将图像数据分解成基向量的集合，且每一井

投影到其上。

@param n_wells –图像分块中的井数。

@param n_frame –图像分块中的井数。

@param image – 与Rawimage结构一样的次序。个别井

帧中的数据，row、col主次序，因此value(row_i,col_j,frame_k) =image[row_i * ncol + col_j + (nrow * ncol * frame_k))

@param n_ sample_ well s –分块的样本中的井数

@param image_sample –向量上方的图像的样本，且如此

@param compressed –输出有损失的经压缩分块

*/

void LossyCompress(float *image);

/***

从经压缩的频域向量重建构图像

数据。

@param n_wells –图像分块中的井数目

@param n_frame –图像分块中帧数。

@param image –与Rawimage结构一样的次序。个别井

帧中的数据，row、col主次序，因此value(row_i, col_j,frame_k) image [row_i * ncol + col_ j + (nrow *neal * frame_k))

*/

void LossyUncompress(float *image);

/***

给定命令行输入，传回基向量的数目

以存储于文件中。

@param _param – 传递到命令行的Dfc参数

*/

static int GetNumCoeff(int _param)

{

return _param;

}

/***

配置平移平滑窗参数。

*/

int SetWindowLength(int length);

inline int GetWindowLength() { return winParamLen; }

float SetWindowAlpha(float alpha);

inline float GetWindowAlpha() { return winParamAlpha; }

float *scratch_pad; //OFT部分结果所需要的中间存储器，当前在LossyCompress()中分配/发布

protected:

void GenerateWindow();

/***

计算产生重要向量所需要的频域相关统计，

其又用以充填

bitsPerFreq向量。

*/

void Emphasis();

ParallelDFT DFT; //离散傅立叶变换处理

object

int n_frames; ///<表示的帧数

int winParamLen;

float winParamAlpha;

std::vector<float> windowCompress; //波纹减少窗

std::vector<float> windowExpand; //在重建构期间的窗使用的逆

};

#endif // DFCCOMPR_H

#include “ParallelDFT.h”

#incl ude <stdio.h>

#incl ude <math.h>

#ifndef PI

#define PI (3.141592653589793f)

#endif // PI

ParallelDFT::ParallelDFT(unsigned int n_pts)

{

SetPts(n_pts);

}

ParallelDFT::~ParallelDFT()

{

m_twiddlesI.clear();

m_twiddlesQ.clear();

m_twiddlesQneg.clear();

}

int Parallel DFT::SetPts(unsigned int n_pts)

{

if (n_pts != m_pts)

{

//仅在改变大小的情况下更新

m_twiddlesI.resize(n_pts);

m_twiddlesQ.resize(n_pts);

m_twiddlesQneg.resize(n_pts);

m_pts = n_pts;

if (n_pts)

{

float kTwiddleScale = -2.0f*PI/static_cast<float>(n_pts);

float fldxVal;

//产生旋转值

for (unsigned int idx = 0; idx < n_pts ; ++idx)

{

fIdxVal = static_cast<float>(idx)*kTwiddleScale;

m_twiddlesI[idx] = cos(fIdxVal);

m_twiddlesQ[Idx] = sin(fIdxVal);

m_twiddlesQneg[idx] = -m_twiddlesQ[idx];

//比重新计算sin(-fIdxVal)快

}

return static_cast<int>(m_pts);

}

/***

计算纯实数输入信号的部分DFT。

*/

int ParallelDFT::PartialDFT(unsigned int n_freqOffset, unsigned intn_numFreqs, const float* srcData, float* dstDataI, float* dstDataQ)

{

int retVal = 0;

unsigned int freqIdx, timeIdx, maxFreq;

register int idx;

register float tmpI , tmpQ;

if (dstDataI && dstDataQ)

{

maxFreq = n_freqOffset + n_numFreqs;

for (freqIdx = n_freqOffset ; freqIdx < maxFreq;++freqIdx)

{

tmpI = 0.0f;

tmpQ = 0.0f;

//此处是发生某一并行处理之处

// SSE/AVX

for (timeldx = 0; timeldx < m_pts ; ++timeIdx)

{

idx (freqIdx * timeldx) % m_pts;

tmpI += srcData[timeIdx] * m_twiddlesI[idx];

tmpQ += srcData[timeIdx] * m_twiddlesQ(idx] ;

}

dstDataI[retVal] = tmpI;

dstDataQ[retVal] = tmpQ;

++retVal;

}

return retVal;

}

/***

计算复频率要素的部分集合的

逆DFT。此处假设绝不会将比

频谱的一半还多的部分回馈到IDFT，这允许利用

傅立叶变换的对称性，同时期望实数结果，避免

大量计算。

*/

int ParallelDFT::PartialIDFT(unsigned int n_freqOffset, unsigned intn_numFreqs, float* dstData, const float* srcDataI, const float* srcDataQ)

{

int retVal = 0;

unsigned int freqIdx, timeIdx, maxFreq;

register int idx;

register float tmp, scale;

if (dstData && n_ numFreqs)

{

maxFreq = n_freqOffset + n_numFreqs;

scale = 1.0f / static_cast<float>(n_numFreqs);

// (1/N)*(N/n_numFreqs)

for (timeldx = 0; timeldx < m_pts; ++timeldx)

{

tmp = 0.0;

//此处发生一些并行处理

// SSE/AVX

for (freqldx = n_freqOffset; freqldx <maxFreq; ++freqldx)

{

idx = (freqldx * timeldx) % m_pts;

//假设结果将为纯实数（复

// 部分将消除），因此仅计算实数分量

tmp += srcDataI\[freqIdx] * m_twiddlesI\[idx] - srcDataQ\[freqIdx] * m_twiddlesQneg\[idx];

}

dstData[retVal] = scale * tmp ;

++retVal;

}

return retVal;

}

#ifndef PARALLELDFT_H

#define PARALLELDFT_H

#include <vector>

/***

离散傅立叶变换实施方案经优化以仅处理

实际上将由DFC算法需要的频率要素。

初始实施方案是不可选的，而是含有注释，其中

可替换SSE/AVX指令以得到性能的提升。

*/

class ParallelDFT

{

public :

ParallelDFT() { m_pts = 0 ; m_twiddlesI.clear(); m_twiddlesQ.clear();

m_twiddlesQneg.clear(); }

ParallelDFT(unsigned int n_pts);

~ParallelDFT();

unsigned int GetPts() { return m_pts; }

int SetPts(unsigned int n_pts);

/***

计算纯实数输入信号的部分DFT。

*/

int PartialDFT(unsigned int n_freqOffset, unsigned int n_numFreqs, const float* srcData, float* dstDataI, float* dstDataQ);

/***

计算复频率要素的部分集合的

逆DFT。

*/

int PartialiDFT (unsigned int n_freqOffset, unsigned intn_numFreqs, float* dstData, const float* srcDataI, const float* srcDataQ);

private:

unsigned int m_pts;

std::vector<float> m_twiddlesI;

std::vector<float> m_twiddlesQ;

std::vector<float> m_twiddlesQneg;

};

#endif //PARALLELDFT_H

关于本发明的测序方面可包括Rothberg等人的第7,948,015号美国专利和Rothberg等人的第2010/0137143号、第2009/0026082号和第2010/0282617号美国专利申请公开案中所描述的一个或多个特征，所述专利皆以全文引用的方式并入本文中。

关于本发明的数据分析方面（例如，处理测量结果，产生经预测信号和使用被称为基底的相位模型建模残差等）可包括Davey等人的第2012/0109598号美国专利申请公开案和Sikora等人的第2013/0060482和2013/0090860号美国专利申请公开案中所描述的一个或多个特征，所述公开案皆以全文引用的方式并入本文中。

在各种实施例中，前述方法的一个或多个方面可至少部分使用现场可编程门阵列（“FPGA”）技术和/或图形处理单元（“GPU”）实施。以下文件皆以全文引用的方式并入本文中：Woods, R.等人的“信号处理系统的基于FPGA的实施方案（FPGA-based Implementation of Signal Processing Systems）”，John Wiley & Sons（2008年）；Gallagher, S.的“将 DSP算法映射到FPGA（Mapping DSP Algorithms Into FPGAs）”，其可在如下地址获得：http://www.ieee.li/pdf/viewgraphs/mapping_dsp_algorithms_into_fpgas.pdf；以及Bartholomä, R.等人的“实施关于FPGA的信号处理算法（Implementing Signal Processing Algorithms on FPGAs）”（德国Pforzheim应用科学大学），其可在如下地址获得：

。

图22为根据本发明的示范性实施例的可经配置为用于执行上文所描述的方法的计算机、系统和/或服务器的计算机的简化功能框图。具体地说，在一个实施例中，如图22中所示，实施上文描述的公开内容的任何计算机、系统和/或服务器可为硬件2200的组合件，包含例如用于包数据通信的数据通信接口2260。平台还可包含用于执行程序指令的中央处理单元（“CPU”）2220，其呈一个或多个处理器的形式。平台通常包含内部通信总线2210、程序存储器和用于待由例如ROM 2230和RAM 2240等平台处理和/或传达的各种数据文件的数据存储器，但系统2200通常经由网络通信2270接收编程和数据。服务器2200也可包含输入和输出端口2250以与例如键盘、鼠标、触摸屏、监视器、显示器等输入和输出装置连接。当然，各种服务器功能可以分配方式实施于数个类似平台上，以分配处理负载。或者，服务器可由一个计算机硬件平台的适当编程来实施。

技术的程序方面可被认为是“产品”或“制品”，其通常呈携载在或体现在一类型的机器可读媒体上的可执行代码和/或相关联的数据的形式。“存储”型媒体包含计算机、处理器等的任何或所有有形存储器或其相关联的模块（例如各种半导体存储器、磁带驱动器、磁盘驱动器等），其可以在软件程序化的任何时间提供非暂时性存储。所有或部分软件有时可以通过因特网或各种其它电信网络通信。举例来说，此通信可使得软件能够从一个计算机或处理器载入到另一计算机或处理器中，例如从移动通信网络的管理服务器或主机载入到服务器的计算机平台中，和/或从服务器载入到移动装置。因此，可以带有软件元件的另一类媒体包括光波、电波和电磁波，例如穿过本地装置之间的物理接口、通过有线并且光学的固定网络以及在各种空中连接上使用。携载所述波的物理元件（例如有线或无线链路、光学链路等）也可以视为带有软件的媒体。如本文中所用，除非限制为非暂时性、有形“存储”媒体，否则例如计算机或机器“可读媒体”的术语是指参与向处理器提供执行指令的任何媒体。

虽然本发明所公开的应用、方法、计算机、服务器、装置和系统被描述为示范性参考计算机应用和发射各种类型的数据，但应了解本发明所公开的实施例可适用于任何环境，例如桌上型或膝上型计算机等。而且，本发明所公开的实施例可适用于等效于HTTP或为HTTP的继承者的任何类型的因特网协议。

从本文中公开的本发明的说明书和实践的考虑，本发明的其它实施例将对所属领域的技术人员显而易见。希望仅将说明书以及实例视为示范性的，其中本发明的真实范围和精神由以下权利要求书来指定。

Claims

1.一种用于压缩测序数据的计算机实施的方法，所述方法包括：

接收与传感器阵列上发生的化学事件相关联的波形数据，所述波形数据包含所述传感器阵列的对应多个位置的多个基于时间的波形；

通过至少一个处理器将所述波形数据的每一基于时间的波形转换成频域谱；

通过所述至少一个处理器基于多个所述频域谱产生关键帧；

通过所述至少一个处理器针对所述频域谱中的每一个计算所述频域谱与所述关键帧之间的差；以及

通过所述至少一个处理器编码所述频域谱与所述关键帧之间的每一经过计算的差；

其中将所述波形数据的每一基于时间的波形转换成频域谱包含：

通过所述至少一个处理器使用积分变换将所述波形数据的每一基于时间的波形变换成频域谱。

2.根据权利要求1所述的方法，其中基于多个频域谱产生所述关键帧包含：

通过所述至少一个处理器平均化所述多个频域谱。

3.根据权利要求1所述的方法，其进一步包括：

通过所述至少一个处理器截断所述多个频域谱，

其中产生所述关键帧是基于经截断的多个频域谱。

4.根据权利要求1所述的方法，其进一步包括：

通过所述至少一个处理器截断所述关键帧；

通过所述至少一个处理器截断所述多个频域谱，

其中计算每一频域谱与所述关键帧之间的所述差包含计算经截断的多个频域谱与经截断的关键帧之间的差。

5.根据权利要求1所述的方法，其进一步包括：

通过所述至少一个处理器基于每一频域谱与所述关键帧之间的经过计算的差来确定编码所述频域谱中的每一个所需要的位数。

6.根据权利要求5所述的方法，其进一步包括：

通过所述至少一个处理器计算所述频域谱与所述关键帧之间的所述经过计算的差中的每一个的缩放向量，所述缩放向量是基于用以编码相应频域谱的所述位数。

7.根据权利要求1所述的方法，其进一步包括：

通过所述至少一个处理器存储所述关键帧和所述频域谱与所述关键帧之间的多个经编码的差。

8.一种用于压缩测序数据的系统，所述系统包括：

数据存储装置，其存储用于压缩测序数据的指令；以及

处理器，其被配置为执行所述指令以执行方法，所述方法包含：

将所述波形数据的每一基于时间的波形转换成频域谱；

基于多个所述频域谱产生关键帧；

针对所述频域谱中的每一个计算所述频域谱与所述关键帧之间的差；以及

编码所述频域谱与所述关键帧之间的每一经过计算的差；

使用积分变换将所述波形数据的每一基于时间的波形变换成频域谱。

9.根据权利要求8所述的系统，其中基于多个频域谱产生所述关键帧包含：

平均化所述多个频域谱。

10.根据权利要求8所述的系统，其中所述处理器被进一步配置为执行所述指令以执行所述方法，所述方法进一步包含：

截断所述多个频域谱，

其中产生所述关键帧是基于经截断的多个频域谱。

11.根据权利要求8所述的系统，其中所述处理器被进一步配置为执行所述指令以执行所述方法，所述方法进一步包含：

截断所述关键帧；

截断所述多个频域谱，

12.根据权利要求8所述的系统，其中所述处理器被进一步配置为执行所述指令以执行所述方法，所述方法进一步包含：

基于每一频域谱与所述关键帧之间的经过计算的差来确定编码所述频域谱中的每一个所需要的位数。

13.根据权利要求12所述的系统，其中所述处理器被进一步配置为执行所述指令以执行所述方法，所述方法进一步包含：

计算所述频域谱与所述关键帧之间的所述经过计算的差中的每一个的缩放向量，所述缩放向量是基于用以编码相应频域谱的所述位数。

14.根据权利要求8所述的系统，其中所述处理器被进一步配置为执行所述指令以执行所述方法，所述方法进一步包含：

存储所述关键帧和所述频域谱与所述关键帧之间的多个经编码的差。

15.一种存储指令的非暂时性计算机可读介质，所述指令在由计算机执行时使得所述计算机执行用于压缩测序数据的方法，所述方法包含：

通过所述至少一个处理器基于多个所述频域谱产生关键帧；