CN113140257A - 一种基因测序信号去除串扰的方法 - Google Patents
一种基因测序信号去除串扰的方法 Download PDFInfo
- Publication number
- CN113140257A CN113140257A CN202010061629.2A CN202010061629A CN113140257A CN 113140257 A CN113140257 A CN 113140257A CN 202010061629 A CN202010061629 A CN 202010061629A CN 113140257 A CN113140257 A CN 113140257A
- Authority
- CN
- China
- Prior art keywords
- pits
- sequencing
- micro
- crosstalk
- signal
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Pending
Links
- 238000012163 sequencing technique Methods 0.000 title claims abstract description 197
- 238000000034 method Methods 0.000 title claims abstract description 86
- 108090000623 proteins and genes Proteins 0.000 title claims abstract description 68
- 238000012216 screening Methods 0.000 claims abstract description 13
- 238000006243 chemical reaction Methods 0.000 claims description 56
- 230000000694 effects Effects 0.000 claims description 9
- 230000002093 peripheral effect Effects 0.000 claims description 5
- 230000008030 elimination Effects 0.000 claims 2
- 238000003379 elimination reaction Methods 0.000 claims 2
- 239000011159 matrix material Substances 0.000 description 24
- 238000004364 calculation method Methods 0.000 description 13
- 238000012545 processing Methods 0.000 description 9
- 238000001514 detection method Methods 0.000 description 8
- 238000012937 correction Methods 0.000 description 6
- 238000000605 extraction Methods 0.000 description 5
- 108020004414 DNA Proteins 0.000 description 4
- 238000010586 diagram Methods 0.000 description 4
- 239000012634 fragment Substances 0.000 description 3
- 238000012165 high-throughput sequencing Methods 0.000 description 3
- 238000012067 mathematical method Methods 0.000 description 3
- 230000005855 radiation Effects 0.000 description 3
- 238000012935 Averaging Methods 0.000 description 2
- 102000001999 Transcription Factor Pit-1 Human genes 0.000 description 2
- 108010040742 Transcription Factor Pit-1 Proteins 0.000 description 2
- 238000002474 experimental method Methods 0.000 description 2
- 108020004707 nucleic acids Proteins 0.000 description 2
- 102000039446 nucleic acids Human genes 0.000 description 2
- 150000007523 nucleic acids Chemical class 0.000 description 2
- 238000012805 post-processing Methods 0.000 description 2
- 102000053602 DNA Human genes 0.000 description 1
- 238000001712 DNA sequencing Methods 0.000 description 1
- 101150054854 POU1F1 gene Proteins 0.000 description 1
- XUIMIQQOPSSXEZ-UHFFFAOYSA-N Silicon Chemical compound [Si] XUIMIQQOPSSXEZ-UHFFFAOYSA-N 0.000 description 1
- 238000002679 ablation Methods 0.000 description 1
- 230000015572 biosynthetic process Effects 0.000 description 1
- 230000000903 blocking effect Effects 0.000 description 1
- 150000001875 compounds Chemical class 0.000 description 1
- 238000012790 confirmation Methods 0.000 description 1
- 238000013075 data extraction Methods 0.000 description 1
- 239000011521 glass Substances 0.000 description 1
- 230000003993 interaction Effects 0.000 description 1
- 238000004020 luminiscence type Methods 0.000 description 1
- 239000000463 material Substances 0.000 description 1
- 230000003287 optical effect Effects 0.000 description 1
- 238000003672 processing method Methods 0.000 description 1
- 230000000191 radiation effect Effects 0.000 description 1
- 239000000376 reactant Substances 0.000 description 1
- 230000001105 regulatory effect Effects 0.000 description 1
- 229910052710 silicon Inorganic materials 0.000 description 1
- 239000010703 silicon Substances 0.000 description 1
- 238000004557 single molecule detection Methods 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
- 238000003786 synthesis reaction Methods 0.000 description 1
- 238000012360 testing method Methods 0.000 description 1
- 238000007671 third-generation sequencing Methods 0.000 description 1
- 235000012431 wafers Nutrition 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
- G16B30/00—ICT specially adapted for sequence analysis involving nucleotides or amino acids
- G16B30/10—Sequence alignment; Homology search
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
- G16B40/00—ICT specially adapted for biostatistics; ICT specially adapted for bioinformatics-related machine learning or data mining, e.g. knowledge discovery or pattern finding
Landscapes
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Health & Medical Sciences (AREA)
- Engineering & Computer Science (AREA)
- Medical Informatics (AREA)
- Biophysics (AREA)
- Theoretical Computer Science (AREA)
- Spectroscopy & Molecular Physics (AREA)
- General Health & Medical Sciences (AREA)
- Evolutionary Biology (AREA)
- Biotechnology (AREA)
- Bioinformatics & Computational Biology (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Data Mining & Analysis (AREA)
- Software Systems (AREA)
- Artificial Intelligence (AREA)
- Evolutionary Computation (AREA)
- Epidemiology (AREA)
- Databases & Information Systems (AREA)
- Public Health (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Bioethics (AREA)
- Chemical & Material Sciences (AREA)
- Analytical Chemistry (AREA)
- Proteomics, Peptides & Aminoacids (AREA)
- Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)
Abstract
本发明提供一种基因测序信号去除串扰的方法。利用孤立数据点进行测算，最终获取复杂环境下的真实测序信号。该方法适用于3端不封闭的测序方法。串扰的方法包含，(1)测序获得测序结果图像；(2)获取测序微坑的信号强度；(3)设定信号强度阈值；(4)筛选出孤立亮坑；(5)估算串扰(6)应用串扰平均值获得校正的测序信号强度值。
Description
技术领域
本发明涉及一种基因测序测序信号去除串扰的方法，属于基因测序领域。
背景技术
基因测序是近年来发展很快的行业。基因测序的过程实际是一个双链DNA分子合成的过程。高通量的基因测序中，由于其信号强度极低，微坑极小，一般在0.2-2微米的范围内进行反应，因此，信号的重复性并不严格。也就是说，在高通量测序的过程中，信号的重复性和一致性并不会很完美。在这个情况下，测序信号的处理是一个很重要的提高测序质量的方式。现代化的信号处理中，可能有很多技术都可以做到这种相对简单的信号。然而，基因测序由于其独特性，并且受到前处理过程的影响比较大，因此，需要开发独特的信号处理方法。基因测序中，微坑或者说数据点，一般都是在0.2-2微米的大小。基因信号检测的过程中，需要用物镜或者其它检测方式获得测序的信号。在多重考虑下，例如使用放大倍率比较大的物镜固然可以获得更好的信号，有更清晰的像素，但是会严重拖慢基因测序的整个流程，一般每个数据点的像素值或者其它信号值并不会很多。在这种情况下，两个数据点之间的相互影响必然会比较严重，后处理的方式也严重影响了基因测序的信号好坏。本发明公开一种基因测序去串扰的方法，利用孤立亮坑计算其对于周围微坑的串扰，从而修正基因测序的信号。
一般的，由于测序信号强度比较低，因此，基因测序的每个数据点的信号强度并不是很高。在这种情况下，背景信号是需要考虑的。测序的信号中掺杂了背景信号也属于本领域的常识。背景信号一般指的是未发生测序反应，也未受到其它信号影响的数据点信号。简单的，当一个测序反应的数据点未受到其它发生测序反应的数据点影响的时候，其实际强度信号可以认为是观测到的强度信号减去背景值。本发明公开一种3端不封闭的基因测序信号去除串扰的方法，利用孤立数据点进行测算，最终获取复杂环境下的真实测序信号。
发明内容
本发明公开一种基因测序信号去串扰的方法，其特征在于包含，
(1)利用基因测序芯片对待测基因序列进行测序反应，获得测序结果图像；
(2)通过测序结果图像获取基因测序芯片上微坑的信号强度；
(3)设定信号强度阈值，高于阈值的信号对应的微坑发生测序反应，低于阈值的信号对应的微坑未发生测序反应；定义发生测序反应的微坑为亮坑，定义未发生测序反应的微坑为暗坑；
(4)筛选出孤立亮坑；
(5)通过孤立亮坑估算亮坑向周围相邻微坑的串扰,取得串扰平均值；
(6)应用串扰平均值获得校正的测序信号强度值；
其中所述的孤立亮坑，指的是该亮坑的直接相邻的周围微坑均为暗坑；
其中步骤(1)中所述的测序反应指的是3端不封闭的基因测序方法。
根据优选的实施方式，与孤立亮坑直接相邻的微坑称为一级串扰微坑，与孤立亮坑隔一个微坑，并且与一级串扰微坑相邻的微坑称为二级串扰微坑；步骤(4)筛选孤立亮坑的时候，选择一级串扰微坑与二级串扰微坑都没有发生反应的孤立亮坑。
根据优选的实施方式，步骤(3)还包括，计算背景信号平均值；通过未发生反应的微坑计算背景信号。
根据优选的实施方式，步骤(4)还包括计算背景信号；通过未发生反应的微坑计算背景信号。
根据优选的实施方式，在去除背景信号影响以后，可以统计串扰平均值。
根据优选的实施方式，所述的统计串扰平均值包括一级串扰微坑的统计串扰平均值，也包括二级串扰微坑的统计串扰平均值。
根据优选的实施方式，所述的测序使用的基因测序芯片的至少一个内表面上有阵列的微坑。
根据优选的实施方式，所述的测序指的是变合成边测序。
根据优选的实施方式，步骤(2)所述通过测序结果图像获取基因测序芯片上微坑的信号强度，指的是参考微坑的物理尺寸设定范围，获取该设定范围内的测序信号强度平均值。
根据优选的实施方式，设定范围以后，处于范围边界上的图像像素点，按照其在设定范围内的比例计算强度。
根据优选的实施方式，所述的参考微坑的物理尺寸设定范围过程中，所设定的范围小于等于微坑的周期。
本发明还公开一种基因测序信号去串扰的方法，其特征在于包含，
(1)利用基因测序芯片对待测基因序列进行测序反应，获得测序结果图像；
(2)通过测序结果图像获取基因测序芯片上微坑的信号强度；
(3)设定信号强度阈值，高于阈值的信号对应的微坑发生测序反应，低于阈值的信号对应的微坑未发生测序反应；定义发生测序反应的微坑为亮坑，定义未发生测序反应的微坑为暗坑；
(4)筛选出孤立暗坑，通过孤立暗坑的信号强度获取背景信号平均值；
(5)筛选出孤立亮坑；
(6)通过孤立亮坑估算亮坑向周围相邻微坑的串扰,并取得串扰平均值；
(7)应用串扰平均值和背景信号平均值去卷积获得校正的测序信号强度值；
其中所述的孤立亮坑，指的是该亮坑的直接相邻的周围微坑均为暗坑；
其中步骤(1)中所述的测序反应指的是3端不封闭的基因测序方法。
本发明公开一种基因测序信号去除串扰的方法，该方法适用于3端不封闭的测序。测序的过程中，每次测序反应可能不止延伸一个碱基。可以出现多倍测序信号。该方法具备以下优势：(1)更加适用于3端不封闭的测序；对于3端封闭测序的效果不明显。(2)通过筛选出来孤立的亮坑，精确获得亮坑对于周围微坑的串扰；(3)相比于不应用去除串扰的方法，测序结果的正确率提升巨大。
除特殊说明外，本发明中所涉及的词均为本领域的通用说法。
附图说明
图1.微坑串扰示意图；
图2.去除串扰流程图；
图3.孤立坑与非孤立坑示意图；
图4.背景坑示意图；
图5.串扰矩阵示意图；
图6.测序图像画圆圈定范围；
图7.信号强度直方图，横轴表示强度，纵轴表示以最大个数归一化后的百分比；
图8.筛选出的亮坑，在亮坑的位置，用“x”的形式进行表示；
图9.孤立亮坑标记图。
具体实施方式
为了进一步说明本发明的核心内容，现将本发明用下面的例子作为说明。实施例是为了进一步解释发明内容部分，并不对于本发明造成限制。
测序是一个对于未知或者已知序列进行检测的过程。基因测序的过程中，一般是一个化学反应的过程。由于基因十分小，单个的基因测序信号是难检测的。目前的单分子检测技术也并不适用于大规模的基因检测。因此，二代测序技术是目前基本通用的测序技术。二代测序技术是将需要检测的基因片段在CPR等过程的基础上，扩增成更多的序列，以便于检测。由于数据量十分大，因此每个数据点的尺寸也是十分微小。
基因测序芯片的空间结构是有规则的，例如四方形或者六方排列等特点的微坑，不能是杂乱随机的微坑。常规的基因测序芯片也都是有规则的。因此，本方法适合于常见的带有微坑的基因测序芯片。带有微坑的基因测序芯片，一般是在测序芯片的一个或者两个底面上，有阵列的、规则的微坑结构。而每个微坑结构都是一个检测点。
在高通量DNA测序中，理想情况下，一个坑内对应一条序列信号强度，但是实际中，由于光学PSF等的影响，如图1所示，每个坑都会影响到周围坑的强度，即每个坑的强度会串扰到周围坑中。在一些常规的测序仪(比如illumina等)中，对DNA信号的检测只是定性的检测，即只需判断信号的有无，在这些仪器中，串扰的影响很小，可不考虑。但是在一些新型测序实现中，需要定量的检测到每个反应坑的信号强度，这时对每个坑真实的反应强度需要较精确的测量，特别是当高DPL信号存在时，串扰将严重影响测序质量。如图1，当坑1中的信号原始DPL为10，坑2中原始信号DPL为1时，若串扰率为10％，则坑2的实际检测信号约为1+10*10％＝2。见图1右侧图片，及时微坑并不直接相邻，也可能产生相互影响的作用。如果不去除这串扰的影响，单个坑提取的强度已不能实际反应单个坑反应所生产的光强，造成最后结果的错误。因此本发明设计一种测序信号中去除串扰的方法，有效去除周围坑串扰的影响。
本发明中涉及到了阈值一词。阈值一词并不是新名词。例如专利CN102834828 B中也提到了阈值的使用。阈值也可以称为临界值。本发明中使用阈值概念，其意义为：高于阈值的信号对应的微坑发生测序反应，低于阈值的信号对应的微坑未发生测序反应。通常的，这个数值是一个统计值或者计算值。阈值的选择并不是复杂的程序。参见附图1，反应的微坑与未反应的微坑的信号有明显差异。那么阈值的选择会比较简单，甚至可以有很多选择。因此，其并不复杂。并且结合之前的技术，阈值是很容易选择的。另外还可以参见实施例1，及其附图7，比较容易的以暗坑和亮坑峰的之间最低值作为阈值，区分是暗坑还是亮坑。这种区分并不会对于结果造成很大的影响。
本发明中所述的基因测序指的是，3端不封闭的基因测序方法。3端封闭的时候，获得的基因测序信号可以表示为1和0。当3端不封闭的时候，由于其测序连续发生，会一次性反应多个碱基，因此，可以获得大于1的信号，例如信号强度为1,2,3,4,5等等。
此处，信号1表示可以简单的认为是测序反应延长一个碱基。在基因测序的过程中，可以认为是一个数据点上面的所有核酸分子，延长了一个碱基。当然，实际反应的时候，可能该数据点上并不是所有的核酸发生反应，这也是测序数据错误的重要原因。
本发明所述的基因测序信号1和0表示的是刨除背景值以后的基因测序信号的相对强度值。例如延长一个碱基的时候，可以认为信号强度为1，延长N个碱基的时候，信号强度为N。实际测序的时候，这种信号强度并不是简单的倍数关系，可能跟很多的影响因素有关。实际测序获得的信号值并不是整数,但是其通过与单碱基测序信号做比较,即可以得到其整数的信号强度。
当信号仅为1和0的时候，本方法的意义不大。因为其只需要判断信号的有和无。当然在刨除背景信号的情况下。
当信号大于1的时候，例如3,6等，由于微坑之间串扰等影响因素，精确判定微坑内的信号值是比较困难的。本发明所述的方法具备很大的优势。
本发明所述的微坑指的是，测序反应中，所使用的微反应室。测序反应的过程中，检测测序信号需要一定强度的测序反应信号。目前的二代测序，都是使用PCR的方法对于基因片段进行扩增，在一定的小范围内，基因片段的数量达到一定程度的时候，可以通过物镜被CCD等检测到。本发明中，所述的微坑指的是，测序反应的过程中，基因片段所在的微小的反应室，一般的反应室的大小在0.2-3微米直径；优选的，在0.3-2微米直径；更优选的，在0.35-1.8微米直径。测序反应室的深度一般的在0.2-3微米；优选的，在0.3-2微米；更优选的，在0.35-1.8微米。微坑属于本领域的常见词语。本发明所述的微坑同本领域的常规含义一致。
本发明中，所述测序指的是高通量的二代基因测序。
本发明所述的微坑，指的是，通过一定的技术手段，在玻璃、硅片等表面形成的，规则排列的微坑。使用微坑属于本领域的常见技术。
基因测序的过程中，所用的测序方法的不同，会影响信号的类型。但是3端不封闭的测序信号，都属于本发明的数据类型范围。
基因测序属于一个复杂的化学反应过程，由于其微坑十分小，信号值十分弱，所以背景信号对于测序信号有很大的影响。简单的，当某个微坑发生反应的时候，测序图像中可以认为是亮坑。亮坑对于周围微坑有信号提升的辐射影响。因此，当一个亮坑是孤立亮坑的时候，其测序真实信号就是测序图像直接提取的信号减去背景信号。背景信号的原因是多方面的，例如使用的基因测序芯片的材料，例如反应物或者其它化合物的残留等等。基因测序的信号提取的过程中，不可避免的受到背景信号的影响。特殊的，3端不封闭的测序反应中，背景信号的影响并不是很严重，例如illumina测序仪所用到的测序方法。3端不封闭的测序信号处理的过程可以比较简单，比如通过经验或者参考序列，或者其它方式，只需要设定一个阈值，则高于该阈值的信号认为是发生测序反应的信号，低于该阈值的信号认为是未发生测序反应的信号。
针对3端封闭的测序信号，其数据类型是典型的有和无。也就是说，可以表示成0和1的形式。这种情况下，本发明所述的方法也是可以用的。比如当数据不够好的时候，比如PCR轮数比较低，或者需要获得更好的数据结果；为了使得测序更加精确，本发明所述的方法也是可以使用的。但是，可以知道的是，这种数据类型比较简单，比如串扰的强度认为是1％到5％，那么数据的提高程度是有限的，并不会实质的提升测序数据的质量。根据简单的实验预测，这种情况下，数据准确性比不使用任何去除串扰的技术可以将错误率从2％降低到1.8％。当然这个数据根据原始的测序错误率还会有所变化。这里必须说明的是，目前所有测序仪公开的准确度数据中，都是经过公开方法的数据计算之后的结果，并不能直接当做原始数据的错误率进行处理。
背景信号的处理方式或者步骤跟其原理有关。本发明中，背景信号处理的步骤是有区别的。例如，获得串扰(信号)平均值之前，可以统计暗坑的信号值，并将该信号值作为背景信号平均值。此处所述的统计，并不要求统计所有的结果，由于测序芯片在小区域内的信号强度比较平行，一般统计部分结果即可。因此，背景信号平均值可以在区分亮暗坑以后获得，或者在筛选出孤立亮坑以后获得。
再例如，当某个微坑属于孤立的暗坑的时候，那么他并不会受到串扰平均值的影响，因此，可以直接统计该区域部分或者全部孤立微坑的信号值，即为背景信号平均值。获取背景信号平均值的过程中，不能仅用一两个微坑的背景信号作为背景信号平均值。区分的，当某个区域的孤立暗坑比较多的时候，统计大部分或者全部的孤立暗坑的信号即可作为背景信号平均值。
针对3端不封闭的测序数据，本发明所述的方法就显得尤为重要，简单的进行对比，当不使用本方法的时候，测序数据的错误率是2％，使用本方法以后，测序数据的错误率可以降低到0.5％，甚至0.1％。此处所列举的数据只是简单的相同实验条件下的结果，并不是最终的数据结果。一般的针对多碱基测序，去除串扰的方法是需要同失相校正等方法联合使用的。根据申请人以前的公开的专利中的实验，完全不进行校正的时候，测序错误率大约为2％；使用一般的简单数据处理，例如常规的数学方法，就可以将数据的错误率大大的降低，例如可以达到99％的正确率。当使用去除串扰的方法的时候，这个数据正确率至少可以达到99.9％，到99.99％，相当于行业Q40标准；当同失相校正方法(201610899880.X，2015109448785)等其它现有技术联合使用的时候，数据正确率可以达到行业Q60的标准。
提取强度的过程是将数据点的测序信号以可信的方式表征的过程。单个数据点的强度提取是一个比较简单的过程。可以利用简单的范围圈定，像素确认即可得到。但是基因测序多倍信号的特点决定了，单纯提出来每个微坑的强度并不能获知微坑的实际或者真实测序信号。范围的划定是一个图像处理问题，一般的，测序微坑是圆形的。画一个圆形即可圈定微坑的范围。当然，为了精确计算，可以将该范围设定成比微坑周期更小的模式。这是根据微坑的尺寸与像素值的关系来操作的。属于一个经验值。
申请人之前公开的技术(201610899880.X)属于对于数据的后处理。配合本发明的数据提取的过程，可以对于测序的数据有更真实的表现。另外，申请人还公开了数据错误校正的方法(cn2015109448785)，后面也称为ECC方法，该方法适合于多碱基的有相互正交关系的测序数据。更进一步的说，本发明所述的方法适合于多碱基测序，201610899880.X所述的方法适合于所有测序领域，ECC方法适合于多碱基的相互正交的测序。因此，当多碱基相互正交测序的时候，本发明所述的方面和前面两种方法是可以联合使用的，能够更进一步的增加数据的准确性。简单的测试结果表明，在多碱基相互正交测序数据处理的时候，不使用本发明所述的方法测序准确度为行业Q40标准，使用本发明所述方法的时候，可以达到行业Q60的标准，200bp。由此也可以知道，本发明所述的方法对于数据的准确性有了较大的改善。
每个微坑的DPL实际是每个微坑的像素值。本发明中，采用了一般的DPL，比如4-64。单个微坑或者数据点的像素值是实际测序需要考察的重要要素。当像素太高的时候，总的数据量会比较低，高通量测序变的没有太多意义。
本发明中，所述的微坑也指的是单个数据点。因此，本发明中，微坑一词也可以更换成数据点。因此，相应的暗坑也是暗数据点，亮坑也是亮数据点。
本发明流程图如图2所示。
首先的需要进行测序然后获得测序的图像。根据前面所述，本发明中的所述测序指的是3端不封闭的测序。
提取强度：
提取强度的过程是对于每个测序点的信号读取的过程。根据微坑数据点的大小，提出出该数据点的信号强度。
对芯片中的每个坑，按照坑在图像上的实际大小，提取出该坑在图像中所反应的实际像素位置，然后在图像中计算这些像素的平均值为该坑的信号强度。这样可以得到每个坑的信号强度。由于芯片上的坑在空间上等间隔阵列排列，因此可以将所有芯片上的坑排列成一个N*M的强度矩阵Iobs，矩阵的每个元素对应每个坑的强度，即强度矩阵的第(i,j)个位置表示的是芯片相应位置的第(i,j)个坑。
计算串扰和背景信号：
提取强度以后，可以明显的发现，部分数据点会明显的比较亮，也有很多数据点明显的比较暗，这些比较亮的数据点所在的微坑就是发生反应的微坑；比较暗的数据点所在的微坑就是没有发生反应的微坑。这个数据拿到以后，就可以用比较常规的方法简单的表示一个数据点是亮还是暗。由于测序信号之间离得比较近，因此，亮的微坑可能会对于周围有个辐射的作用，也就是常常说的信号串扰。因此，当一个比较暗的微坑周围，比如两层微坑的范围内，都没有亮的微坑的时候，可以认为这个暗的微坑属于背景信号。统计很多的背景信号，给出一个平均值，就可以获得平均的背景信号。计算串扰的过程相当于计算热辐射的过程；不同的是，辐射的过程中，一个亮点的周围并不是一个平均衰减的过程，而是以微坑为单位的剧烈衰减。一般的，亮点影响1-2层的周围数据点，更多的范围并不会受到影响。因此，只需要计算一个亮点周围的微坑亮度，并且除去背景的影响，就可以简单的获得串扰的数据。当然在这个筛选的过程中，如果多个亮点相互影响，会导致计算不方便，因此，选择孤立的亮点，也就是说这个亮点周围很大的区域都没有其它的亮点，就很容易的获得亮点对于周围的辐射或者说串扰。
首先筛选亮坑。对芯片中的每个坑，其强度有两种形式，若坑中有DNA片段，则会有测序反应，将会发光，其信号较强；反之，若没有DNA片段，则不会发光，其信号较弱，因此通过对所有强度值做直方图分布，直方图上至少会有两个峰，一个是不发光的暗坑形成的强度峰，一个是发光的亮坑形成的强度峰。在这两个峰之间设定阈值，即可筛选出所有亮坑。
然后在所有亮坑中，挑选孤立坑。如图3所示。孤立坑的定义为该亮坑周围x方向和y方向距离一定范围内都没有其他亮坑。该范围的大小事先根据经验确定。在实现中，对每个亮坑，我们统计所有其它亮坑到其的距离，若这些距离中，最小距离大于事先设定的范围，则标记其为孤立坑。
再次，估计背景信号B。对上一步中的亮度坑，我们标定其附近一定范围内的坑都为干扰坑，如图4中方块标记区域所示。除去所有亮度坑和干扰坑，余下的所有坑我们标记为背景坑。统计这些背景坑的信号平均值，作为背景信号强度B(背景信号平均值)。具体的计算的时候，也可以有多种选择，背景坑的信号平均值并不一定需要所有的数据，比如挑选多个未发生测序反应的点进行统计平均也是可以的。例如，当一个未发生反应的微坑直接相邻的一级相邻微坑，以及与一级相邻微坑相邻的二级相邻微坑都没有发生反应的时候，这个未发生反应的微坑的信号也代表了背景信号；多个背景信号的平均值就代表了背景信号平均值。计算背景信号平均值的实际意义就是统计未发生反应微坑的背景值。
最后，计算串扰矩阵F。计算强度矩阵其实是把信号以矩阵的形式表示出来，便于后面的计算。也可以用其它方式表示。我们通过孤立坑来计算各个方向的串扰值。每个坑会串扰影响周围一系列坑的的信号强度。以一级串扰为例，会影响周围八个坑的强度。对每个方向的串扰值，其计算方法为：(该坑强度-背景信号)/(亮坑强度-背景信号)。例如有强度矩阵
其中a22强度为筛选出的孤立亮坑，B为背景信号，则串扰矩阵为：
最后对所有孤立坑计算得到的串扰矩阵去平均，可得到最后的串扰矩阵为
其数据方位参见图5。
校正信号：
校正信号的过程就是去除串扰和背景影响的过程。单个的信号计算比较简单。大量数据的时候，可以使用去卷积的方式很快得到校正的信号。
假设Iori为真实信号(去串扰后信号)，则信号模型可表示为
其中
表示图像卷积操作,B为背景信号，Iobs为检测得到的信号，F为串扰矩阵。故真实信号为上述方差的解，即
Iori＝deconv(Iobs-B，F)
其中，deconv(.)表示解卷积操作。
解卷积在声学、图像处理等领域中广泛应用。
必须说明的是，本发明中所有的计算方式或者公式均只是遵循了测序原理的基础概念。也就是说，所有的公式并没有带来额外的意义，仅仅是测序数据的表示方式。这种表示方式可以是多种多样的，例如简单的一个数据点，可以用手算的方式根据其本发明指出的方法进行计算。因此，公式只是本发明方法的简单替代而已。可以有其它方法进行简单的计算。
本发明中，所述的测序指的是二代测序，并不是三代测序。
本发明中，所述的信号强度指的是一般意义上图像上提取的强度信息。属于本领域的常规用语。
本发明中，串扰一词属于其基本含义。在测序上指的是一个微坑或者数据点的信号对于周围微坑或数据点的影响。
本发明中，已经定义了一级串扰微坑。简单的，一个数据点的周围的其余直接相邻的8个数据点属于其一级串扰微坑或者一级串扰数据点。数据点从图像上直观的来说是一个二维点的类型。
本发明中定义了二级串扰微坑。简单的，同一级串扰微坑相邻的，第二层的微坑属于二级串扰微坑，一共有16个。
一般的测序的时候，微坑的大小大约为0.5-10微米的量级。本发明所述的测序中，所涉及的微坑为0.5-4微米，优选1-3微米直径，更优选1.5-2.5微米直径。
一般的测序的时候，每个微坑拍照的时候，像素个数需要超过4。本发明所述的测序中，每个微坑的像素个数为3-100，优选12-75，更优选24-64。部分本发明的实施例中用了48个像素点，更低一级或者高一级并不影响本发明所述的方法。
实施例1
在本实施例中，坑位置在二维平面上等间隔分布，坑间距DL＝8个像素。
1.提取强度
根据坑位置和尺寸，如图6所示，画出一系列圆(也可以是矩形等形状)，计算每个圆内所有像素的平均值作为该坑信号强度。由于坑在芯片上横平竖直分布，所以根据坑的几何位置分布，可以得到强度矩阵I，矩阵的每一个元素值表示某一个坑的强度值。
计算信号强度的时候，由于一般的微坑是圆形的，因此采用画圆的方式圈定微坑的位置。当微坑是其它形状的时候，也可以改变该方法。
四个相邻微坑计算信号强度的时候，其共有的相邻非圆形区域可以认为是非计算区域。实际计算的时候，该区域的影响可以忽略，不进行计算。
当选择的圆内像素点非整数的时候，只需要简单的切割即可。例如31.4个像素值，只需要将边缘的一个像素值乘0.4即可简单操作。
2.筛选亮坑
筛选亮坑即为将反应的和未反应的微坑区分开来。利用强度矩阵，一般的数学方法可以简单的区分反应和未反应的微坑。将亮坑和暗坑区分开，也是一个去背景的过程。当一个微坑中未发生化学反应，并且其一级串扰微坑，二级串扰微坑都未发生反应的时候，属于背景微坑。
对所有微坑强度做直方图，如图7所示，以暗坑和亮坑峰的之间最低值作为阈值，区分是暗坑还是亮坑。该阈值是一个信号强度的值，高于阈值则为亮坑，低于阈值则为暗坑。图8为筛选出的亮坑，在亮坑的位置，用“x”的形式进行表示。
3.确定孤立亮坑
孤立亮坑指的是定义该亮坑的周围微坑，以及其二级微坑都未发生测序反应的微坑为孤立亮坑。
简单的说，该亮坑周围的一级相邻微坑以及二级相邻微坑都没有发生反应，则定义该亮坑为孤立亮坑。
对每个亮坑点，定义其坑距为其他所有亮坑点到该点的距离的最小值
其中disti为第i个亮坑的坑距，doti为其位置，|·|表示两个点的距离。
在所有的dist中，如果dist>4*DL，则标记其为孤立坑。
如图9所示，x为筛选出的孤立坑，标记o的为暗坑
4.估计串扰矩阵
对每个孤立亮坑，从1中计算的强度中提取以其为中心的5*5的强度矩阵，上图中对应孤立亮坑的强度矩阵为
| 158.9086 | 155.112 | 154.9645 | 154.5978 | 154.4302 |
| 158.8925 | 164.9705 | 170.1037 | 158.4598 | 155.6899 |
| 157.3772 | 176.9054 | 1219.421 | 179.7076 | 160.4912 |
| 163.0005 | 158.0009 | 181.5133 | 163.2679 | 161.3172 |
| 157.6925 | 158.0391 | 159.7521 | 154.154 | 160.4775 |
对每个强度矩阵，最外圈的16个点为其背景坑，取这16个点的中位数为背景坑亮度，上例中，背景坑亮度B＝157.8658。该处的背景坑信号值得是背景信号平均值。
实际的测序数据中，亮坑的信号强度明显高，但是其直接相邻微坑，也就是一级串扰微坑与其它非相邻微坑相比，信号强度并不会特别大。但是串扰是存在的。因此，串扰需要选择多个孤立亮坑进行估计。
3*3串扰矩阵的每个点计算方式为：(该坑强度-背景信号)/(亮坑强度-背景信号)，这里背景信号B＝157.8658，亮坑强度＝1219.421，因此可以计算得到最终的串扰矩阵：
将所有孤立亮坑的平均串扰矩阵作为最后的串扰矩阵输出F。其中，F矩阵中亮坑强度为1，周围坑的数值为串扰的百分比。
5.去串扰
采用解卷积实现。解卷积或者去卷积是类似领域常见的方法。
例如在MATLAB实现中，可用deconvreg函数实现去卷积效果
J＝deconvreg(I,F)
其中I为原始信号强度，F为前面估计的串扰矩阵，J为去除串扰信号后的强度。经过计算,可以获得所有亮点的强度信号。
基因测序的数据量非常大，本发明仅列举其中部分小区域的处理过程作为例子。
申请人同时使用了专利CN201510822361.9所述的基因测序方法，结合了专利CN201510944878.5以及CN201610899880.X及常规手段，可以将测序数据的准确性提升到行业Q60的标准。
根据申请人以前的公开的专利中的实验，完全不进行校正的时候，测序错误率大约为2％；使用一般的简单数据处理，例如常规的数学方法，illumina等公司已经公开，就可以将数据的错误率大大的降低，例如可以达到99％的正确率。当使用去除串扰的方法的时候，这个数据正确率至少可以达到99.9％，到99.99％，相当于行业Q40标准；当同失相校正方法(201610899880.X，2015109448785)等其它现有技术联合使用的时候，数据正确率可以达到行业Q60的标准。
以上所述，仅为本发明较佳的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，可轻易想到的变化或替换，都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应以所述权利要求的保护范围为准。
Claims (10)
1.一种基因测序信号去串扰的方法，其特征在于包含，
(1)利用基因测序芯片对待测基因序列进行测序反应，获得测序结果图像；
(2)通过测序结果图像获取基因测序芯片上微坑的信号强度；
(3)设定信号强度阈值，高于阈值的信号对应的微坑发生测序反应，低于阈值的信号对应的微坑未发生测序反应；定义发生测序反应的微坑为亮坑，定义未发生测序反应的微坑为暗坑；
(4)筛选出孤立亮坑；
(5)通过孤立亮坑估算亮坑向周围相邻微坑的串扰，取得串扰平均值；
(6)应用串扰平均值获得校正的测序信号强度值；
其中所述的孤立亮坑，指的是该亮坑的直接相邻的周围微坑均为暗坑；
其中步骤(1)中所述的测序反应指的是3端不封闭的基因测序方法。
2.根据权利要求1所述的方法，其特征在于，与孤立亮坑直接相邻的微坑称为一级串扰微坑，与孤立亮坑隔一个微坑，并且与一级串扰微坑相邻的微坑称为二级串扰微坑；步骤(4)筛选孤立亮坑的时候，选择一级串扰微坑与二级串扰微坑都没有发生反应的孤立亮坑。
3.根据权利要求1-2任一项所述的方法，其特征在于，步骤(3)还包括，计算背景信号平均值；通过未发生反应的微坑计算背景信号。
4.根据权利要求1-2任一项所述的方法，其特征在于，步骤(4)还包括计算背景信号；通过未发生反应的微坑计算背景信号。
5.根据权利要求1-4任一项所述的方法，其特征在于，在去除背景信号影响以后，可以统计串扰平均值。
6.根据权利要求5所述的方法，其特征在于，所述的统计串扰平均值包括一级串扰微坑的统计串扰平均值，也包括二级串扰微坑的统计串扰平均值。
7.根据权利要求1所述的方法，其特征在于，步骤(2)所述通过测序结果图像获取基因测序芯片上微坑的信号强度，指的是参考微坑的物理尺寸设定范围，获取该设定范围内的测序信号强度平均值。
8.根据权利要求7所述的方法，其特征在于，设定范围以后，处于范围边界上的图像像素点，按照其在设定范围内的比例计算强度。
9.根据权利要求7所述的方法，其特征在于，所述的参考微坑的物理尺寸设定范围过程中，所设定的范围小于等于微坑的周期。
10.一种基因测序信号去串扰的方法，其特征在于包含，
(1)利用基因测序芯片对待测基因序列进行测序反应，获得测序结果图像；
(2)通过测序结果图像获取基因测序芯片上微坑的信号强度；
(3)设定信号强度阈值，高于阈值的信号对应的微坑发生测序反应，低于阈值的信号对应的微坑未发生测序反应；定义发生测序反应的微坑为亮坑，定义未发生测序反应的微坑为暗坑；
(4)筛选出孤立暗坑，通过孤立暗坑的信号强度获取背景信号平均值；
(5)筛选出孤立亮坑；
(6)通过孤立亮坑估算亮坑向周围相邻微坑的串扰,并取得串扰平均值；
(7)应用串扰平均值和背景信号平均值去卷积获得校正的测序信号强度值；其中所述的孤立亮坑，指的是该亮坑的直接相邻的周围微坑均为暗坑；其中步骤(1)中所述的测序反应指的是3端不封闭的基因测序方法。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010061629.2A CN113140257A (zh) | 2020-01-20 | 2020-01-20 | 一种基因测序信号去除串扰的方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010061629.2A CN113140257A (zh) | 2020-01-20 | 2020-01-20 | 一种基因测序信号去除串扰的方法 |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN113140257A true CN113140257A (zh) | 2021-07-20 |
Family
ID=76808849
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202010061629.2A Pending CN113140257A (zh) | 2020-01-20 | 2020-01-20 | 一种基因测序信号去除串扰的方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN113140257A (zh) |
Citations (7)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CA2367413A1 (en) * | 1999-04-08 | 2000-10-12 | Sir Mortimer B. Davis Jewish General Hospital | Quantitative assay for expression of genes in microarray |
| JP2005149266A (ja) * | 2003-11-18 | 2005-06-09 | Olympus Corp | 画像処理装置、画像処理方法及び画像処理プログラム |
| CN101799916A (zh) * | 2010-03-16 | 2010-08-11 | 刘国传 | 基于贝叶斯估计的生物芯片图像小波去噪方法 |
| US20110134288A1 (en) * | 2009-06-30 | 2011-06-09 | Masanori Kasai | Image processing device and image processing method, imaging apparatus, and computer program |
| CN102354398A (zh) * | 2011-09-22 | 2012-02-15 | 苏州大学 | 基于密度中心与自适应的基因芯片处理方法 |
| CN108070525A (zh) * | 2015-11-19 | 2018-05-25 | 赛纳生物科技(北京)有限公司 | 基因测序芯片 |
| CN110184325A (zh) * | 2018-02-22 | 2019-08-30 | 张家港万众一芯生物科技有限公司 | 基于微孔阵列芯片的单分子文库pcr扩增的基因测序方法 |
-
2020
- 2020-01-20 CN CN202010061629.2A patent/CN113140257A/zh active Pending
Patent Citations (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CA2367413A1 (en) * | 1999-04-08 | 2000-10-12 | Sir Mortimer B. Davis Jewish General Hospital | Quantitative assay for expression of genes in microarray |
| JP2005149266A (ja) * | 2003-11-18 | 2005-06-09 | Olympus Corp | 画像処理装置、画像処理方法及び画像処理プログラム |
| US20110134288A1 (en) * | 2009-06-30 | 2011-06-09 | Masanori Kasai | Image processing device and image processing method, imaging apparatus, and computer program |
| CN102124744A (zh) * | 2009-06-30 | 2011-07-13 | 索尼公司 | 图像处理装置、图像处理方法、成像设备及计算机程序 |
| CN101799916A (zh) * | 2010-03-16 | 2010-08-11 | 刘国传 | 基于贝叶斯估计的生物芯片图像小波去噪方法 |
| CN102354398A (zh) * | 2011-09-22 | 2012-02-15 | 苏州大学 | 基于密度中心与自适应的基因芯片处理方法 |
| CN108070525A (zh) * | 2015-11-19 | 2018-05-25 | 赛纳生物科技(北京)有限公司 | 基因测序芯片 |
| CN110184325A (zh) * | 2018-02-22 | 2019-08-30 | 张家港万众一芯生物科技有限公司 | 基于微孔阵列芯片的单分子文库pcr扩增的基因测序方法 |
Non-Patent Citations (1)
| Title |
|---|
| 张瑜: "基因芯片图像的处理和分析方法研究", 《红外与激光工程》, vol. 35, pages 219 - 222 * |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| US7302348B2 (en) | Method and system for quantifying and removing spatial-intensity trends in microarray data | |
| US20200250851A1 (en) | Image analysis useful for patterned objects | |
| WO2018068600A1 (zh) | 图像处理方法及系统 | |
| US10403000B2 (en) | Methods and systems for analyzing biological reaction systems | |
| US7330606B2 (en) | Method and system for extracting data from surface array deposited features | |
| CN110490836B (zh) | dPCR微阵列图像信息处理方法 | |
| EP3843034A1 (en) | Method and device for detecting bright spots on image, and computer program product | |
| EP1450304A1 (en) | Image processing apparatus and method | |
| US7136517B2 (en) | Image analysis process for measuring the signal on biochips | |
| CN109872308B (zh) | 一种矫正微滴式数字pcr通道间微滴位置的方法 | |
| CN114549600A (zh) | 一种荧光图像配准方法 | |
| US20180315187A1 (en) | Methods and systems for background subtraction in an image | |
| US7099502B2 (en) | System and method for automatically processing microarrays | |
| CN113140257A (zh) | 一种基因测序信号去除串扰的方法 | |
| CN117252786A (zh) | 一种基因检测数据增强处理方法 | |
| EP2728502A1 (en) | Method and computer program product for genotype classification | |
| CN114300047A (zh) | 一种基因测序信号强度获取的方法 | |
| EP3254081B1 (en) | Methods and systems for determining optical regions of interest | |
| CN113670865B (zh) | 分辨率板、分辨率评估方法及相关设备 | |
| WO2009126495A2 (en) | Method and system for processing microarray images | |
| Deepa et al. | Automatic segmentation of DNA microarray images using an improved seeded region growing method | |
| Novikov et al. | A robust algorithm for ratio estimation in two-color microarray experiments | |
| US20240177807A1 (en) | Cluster segmentation and conditional base calling | |
| Arteaga-Salas | 9 Image Processing of Affymetrix Microarrays | |
| Zhao et al. | Microarray images processing based on mathematical morphology |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination |