本发明属于生物技术领域，具体涉及一种高通量测序中假阳性变异位点的筛选方法。
背景技术：
随着二代高通量测序技术的普及，近几年获取了大量的基因组重测序数据，其中包括了肿瘤相关的cosmic(catalogueofsomaticmutationsincancer，cosmic)数据库和千人基因组数据库(1000genomeproject)，已经积累了上万个重测序数据样本信息。虽然当前有许多软件或分析流程可以用来对二代测序数据进行snv分析，但是测序和流程分析中存在的许多问题还亟待解决，如检测获得的snv的准确性。通过优化分析流程，设置过滤和比对的参数信息，也不能确定获得的所有snv是正确的。在研究中发现，已有的snv数据库中存在一些可能是错误的突变位点信息，这些潜在的错误snv位点对后续的数据分析和数据挖掘可能会造成影响，导致研究结论可信度降低。
由于基因组序列的复杂性程度比较高，导致我们得到的snv信息可能存在假阳性。例如，人类基因组上的假基因序列和相对应的蛋白质编码基因之间存在高度的相似性区域，其中一些序列不属于高度重复序列，但它们的一些核苷酸序列比较相似。由于当前高通量测序检测到的数据长度为100bp～150bp，如果两段相似性核苷酸序列只有一个核苷酸序列的差异，而这两段序列的长度大于二代测序的最大长度，那么我们很难判断该位点是否发生了变异，该位点的测序分析结果存在两种可能性：(1)存在snv位点，判定该位点为变异位点，则该变异位点的信息可能为假阳性；(2)不存在snv位点，判定该位点没有发生变异，则该位点真实的变异信息被忽略。这种情况使得分析工具，纵使设置严格的参数来减少可能存在的错误，甚至当我们标记重复序列，在分析snv时过去掉重复序列区域的信息，仍然不能避免这类错误的发生。
技术实现要素：
本发明的目的是提供一种高通量测序中假阳性变异位点的筛选方法，能够准确的筛选出高通量测序中由基因组高度相似性序列造成的假阳性变异位点。
本发明所采用的技术方案是：一种高通量测序中假阳性变异位点的筛选方法，按照以下步骤具体进行：
步骤1、下载参考基因组序列，通过计算机语言脚本对上述参考基因组序列的正链信息进行覆盖，覆盖时每次向前步移1bp，获取模拟的二代测序数据，得到了100bp的二代测序reads；
步骤2、将步骤1中获得的reads比对到参考基因组序列上，从比对的结果中查找出99个核苷酸全部匹配、且只有一个错配信息的位点信息，这些位点即为初始的潜在假阳性snv位点；
步骤3、获得步骤2中假阳性snv位点的对应的反链信息，正链信息和反链信息合并得到潜在的假阳性snv位点数据库；
步骤4、将步骤3中获得的数据库位点与cosmic数据库中所有单个替换突变的位点信息进行比较，比较后将相同的位点信息保留，即获得在cosmic数据库中潜在的假阳性变异位点。
本发明的特点还在于，
步骤1中获取模拟的二代测序数据后，去除数据中含有n的序列。
步骤2中比对采用blat比对，所述blat比对结果包括含有gap、错配和完全匹配的信息。
将步骤2中频次低于10的潜在假阳性snv位点过滤掉，再将潜在假阳性snv位点中重复区域位点过滤掉。
对步骤4得到的假阳性位点通过以下方法进行验证：对于步骤4假阳性位点中选择研究目标一致的基因，设计特异性引物进行sanger测序进行验证，即可以筛选出正确的高通量测序中假阳性位点。
步骤1中的参考基因组序列为人的hg19序列。
本发明的有益效果是：本发明首先通过数据模拟人类参考基因组hg19，获得覆盖所有在二代测序中可能检测到的reads，然后将这些reads通过与hg19参考基因组进行比较，获得只有一个位点错配的比对信息，通过层层过滤，最终获得在测序数据中可能为假阳性位点的snv位点信息，将这些位点与cosmic数据集进行比较，发现有6,078个位点实际发生了高频率的假阳性突变，可能会影响研究结果。这些候选错误突变位点希望引起人类基因组研究人员的注意，并为二代测序研究提供候选错误位点，以便更好地获得更准确的结论。
具体实施方式
下面将结合具体实施方式对本发明作进一步阐述。
本发明一种高通量测序中假阳性变异位点的筛选方法，具体按照以下方法进行：
步骤1、模拟数据的准备：
步骤1.1、下载人的hg19参考基因组序列，每一条序列长度设置为100bp，通过perl脚本或python脚本对上述hg19参考基因组序列进行瓦片式覆盖，每次向前步移1bp，获取模拟二代测序数据，即得到长度为100bp的核苷酸序列。
步骤1.2、hg19参考基因组本身含有一些未知的序列区域使用n来替代，在模拟数据中也出现了部分序列含有n序列，因此将步骤1.1中获得的含有一个或者多个n的核苷酸序列给予去除。
步骤2、模拟数据的比对和错配位点信息的获取：
步骤2.1、使用序列比对软件blat将步骤1获得的模拟二代测序数据比对到参考基因组hg19上。
步骤2.2、由于模拟的参考序列都源自参考基因组hg19，因此理论上所有模拟的核苷酸序列和hg19基因组区域完全匹配，blat结果包含了所有可能匹配的信息，除了完全匹配信息，还有一部分比对结果是包含gap，一个或多个错配等信息，本申请分析的目的是找到100bp的序列比对中，只有1个错误匹配，99个完全匹配的序列位点，因此对于步骤2.1中获得的psl文件结果，筛选出只有1个错误匹配且其他位点完全匹配的序列，这些位点即为初始的潜在假阳性snv位点。
步骤2.3、通过位点信息和错位匹配位置信息，计算获得步骤2.2中错误匹配位点在基因组上的位置信息，共获得7,917,212,459个位点，其中有3,965,877,977个位点排布的负链上，有3,951,334,482个位点在正链上。
步骤3、潜在假阳性位点数据库的构建：
步骤3.1、重测序中常用30x的测序，设置的reads数值为3，即至少有3条序列覆盖snv位点，对于该测试中模拟的数据为100x的数据量，因此设置过滤频次为>＝10则认为是在现实测序中可能被检测到的snv位点，低于10则认为是偶然性获得的snv位点。
步骤3.2、在实际的分析过程中，如果snv位于重复序列区域，则会选择屏蔽掉。该方法模拟获取的位点，为基因组上相似区域的序列，其中包含了基因组重复区域的序列，因此将步骤3.1中获得的位点信息中，位于重复区域的去除。人基因组中重复区域的位点占步骤3.1中获得位点的98.35％。
步骤3.3、一开始为了减少数据计算量的问题，本发明只模拟了hg19基因组的正链信息，没有模拟互补链信息，为了保证获得潜在假阳性位点数据库的完整性，对步骤3.2中获得的位点信息进行互补配对并合并，最终获得了1,142,718个位点作为最终人基因组中潜在假阳性位点数据库。
步骤4、cosmic数据库中潜在假阳性变异位点:
步骤4.1、下载cosmic数据信息，有5,962,444个突变信息，本发明提取了所有的单个替换突变，其中包括1,302,143个同义突变，3,898,856个错义突变和258,529个无义位点信息。最后，一共获得了2,915,899个突变位点以备后续分析。
步骤4.2、将步骤4.1获得的位点变异信息与步骤3.3获得的位点信息进行比较，比较后将相同的位点信息保留，发现步骤4.1中共有4,695个位点出现在潜在的假阳性变异位点数据库中。
步骤4.3、对于步骤4.2中获得的变异位点在cosmic中出现的频次进行统计，发现pik3ca(chr3:178936092，c.1634a>c)在cosmic中出现了232次，对变异位点所在的基因组上的核苷酸序列进一步的研究发现，该区域上的序列chr3:178935924-178936114和chr22:17052844-17053034只存在chr3:178936092上a>c的差异。两段相似区域只有一个碱基的差异，导致了在分析过程中该位点被认为是变异位点，该变异位点的真实性需要进一步的测序验证。
步骤5、千人基因组中潜在的假阳性变异位点和实验验证：
步骤5.1、下载千人基因组数据信息，同cosmic数据一样，保留只有单个位点的变异信息，共获得81,653,252个位点信息。
步骤5.2、将步骤5.1中的位点信息与步骤3.3中获得的位点信息进行比较，发现千人基因组数据中共有41,415个位点出现在步骤3.3获得的位点中。
步骤5.3、步骤5.2中的位点信息在千人基因组中的排布特征，发现这些位点可能为假阳性变异位点，例如，rs13329475是一个snp位点位于基因组chr15:25537567(t>c)上，该位点所在的区域与chr15:25537446-25537624具有类似的区域chr8:65202793-65202971，且只有chr15:25537567位点一个碱基的差异。
步骤5.4、通过在chr15:25537567两端设计特异性引物，对亚洲群体的样本，进行pcr扩增并进行一代测序，结果证实该位点并没有发生变异。
本发明的筛选方法中模拟二代测序数据的长度为100bp，并非否定选择150bp的长度，要根据实际测序的长度进行相应的设置，之所以没有选择150bp的读码长度，主要是希望提高覆盖率，从而获得尽可能全面的假阳性snv文库，以便于同cosmic和千人基因组数据进行比较。
对于高度相似区域的序列，在研究变异时会由于参数的设置漏掉一些重要的snp位点信息或得到一些可能错误的snp位点。因此，在进行下一步研究之前，可以利用该发明方法首先获得潜在的错误snv数据库，对获得的snv结果进行比对，并通过sanger测序战略来验证，这可以帮助我们提高分析结果的准确性。另一方面，如果在研究中发现了一些相关的表型和可能导致这种表型的潜在基因，但是在高通量测序中没有找到该基因的变异位点信息，可以通过建立潜在的变异位点库，检索目的基因是否在数据库中出现，如果存在则可以进一步使用sanger测序策略在样本中进行验证。总之，本发明方法可以提供一个非常全面的潜在假阳性位点变异数据库，不仅可以用来过滤二代测序获得的snv结果，提高精确度，而且提醒研究者更加重视这个问题，以获得更准确的结论。
目前现有的筛选过滤假阳性变异位点的方法主要为，1：设置高阈值的参数，阈值设置超过一定的数值，则会导致过滤掉的数据增多，有效数据量减少，其中也会过滤掉一部分阳性数据；2：多种分析流程获得共有的变异位点，在使用两种分析流程计算时，直接增加了一倍的工作量，即流程使用的越多，计算量越大。尤其在样本数量比较多的时候，工作量会非常大。本发明的筛选方法聚焦基因组本身的复杂序列特征，通过数据模拟获得全覆盖的高通量测序数据，通过序列比对和层层过滤分析，找到只有一个碱基差异的相似性序列，最终获得由于该类型序列导致的潜在的假阳性突变位点库。通过对基因组的全覆盖进行一次数据模拟计算，就可以获得完整的假阳性变异位点数据库。只需要对获得的结果文件中的位点进行比较分析即可。减少了运算量的同时，可以在一定程度上提高数据的准确性。