WO2013097328A1 - 基因组indel位点标记方法和装置 - Google Patents

基因组indel位点标记方法和装置 Download PDF

Info

Publication number
WO2013097328A1
WO2013097328A1 PCT/CN2012/071329 CN2012071329W WO2013097328A1 WO 2013097328 A1 WO2013097328 A1 WO 2013097328A1 CN 2012071329 W CN2012071329 W CN 2012071329W WO 2013097328 A1 WO2013097328 A1 WO 2013097328A1
Authority
WO
WIPO (PCT)
Prior art keywords
sequencing
rad
sequence
tag
sequences
Prior art date
Application number
PCT/CN2012/071329
Other languages
English (en)
French (fr)
Inventor
郑泽群
陶晔
汪健
王俊
杨焕明
Original Assignee
深圳华大基因科技服务有限公司
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by 深圳华大基因科技服务有限公司 filed Critical 深圳华大基因科技服务有限公司
Publication of WO2013097328A1 publication Critical patent/WO2013097328A1/zh

Links

Classifications

    • C—CHEMISTRY; METALLURGY
    • C12—BIOCHEMISTRY; BEER; SPIRITS; WINE; VINEGAR; MICROBIOLOGY; ENZYMOLOGY; MUTATION OR GENETIC ENGINEERING
    • C12Q—MEASURING OR TESTING PROCESSES INVOLVING ENZYMES, NUCLEIC ACIDS OR MICROORGANISMS; COMPOSITIONS OR TEST PAPERS THEREFOR; PROCESSES OF PREPARING SUCH COMPOSITIONS; CONDITION-RESPONSIVE CONTROL IN MICROBIOLOGICAL OR ENZYMOLOGICAL PROCESSES
    • C12Q1/00—Measuring or testing processes involving enzymes, nucleic acids or microorganisms; Compositions therefor; Processes of preparing such compositions
    • C12Q1/68—Measuring or testing processes involving enzymes, nucleic acids or microorganisms; Compositions therefor; Processes of preparing such compositions involving nucleic acids
    • C12Q1/6813—Hybridisation assays
    • C12Q1/6827—Hybridisation assays for detection of mutation or polymorphism
    • C—CHEMISTRY; METALLURGY
    • C12—BIOCHEMISTRY; BEER; SPIRITS; WINE; VINEGAR; MICROBIOLOGY; ENZYMOLOGY; MUTATION OR GENETIC ENGINEERING
    • C12Q—MEASURING OR TESTING PROCESSES INVOLVING ENZYMES, NUCLEIC ACIDS OR MICROORGANISMS; COMPOSITIONS OR TEST PAPERS THEREFOR; PROCESSES OF PREPARING SUCH COMPOSITIONS; CONDITION-RESPONSIVE CONTROL IN MICROBIOLOGICAL OR ENZYMOLOGICAL PROCESSES
    • C12Q1/00—Measuring or testing processes involving enzymes, nucleic acids or microorganisms; Compositions therefor; Processes of preparing such compositions
    • C12Q1/68—Measuring or testing processes involving enzymes, nucleic acids or microorganisms; Compositions therefor; Processes of preparing such compositions involving nucleic acids
    • C12Q1/6809—Methods for determination or identification of nucleic acids involving differential detection

Definitions

  • the invention belongs to the technical field of genetic engineering, and in particular relates to a method and device for marking a genome INDEL (Insert-Deletion) site. Background technique
  • the genomic INDEL site marker refers to the difference in the whole genome between the two parents, and one parent has a certain number of nucleotide insertions or deletions in the genome of the parent relative to the other parent.
  • the acquisition of INDEL locus information can have many important applications, such as building genetic maps, genotyping, molecular marker breeding, disease detection, and more.
  • Second-generation DNA sequencing technology is a high-throughput, low-cost sequencing technology.
  • the basic principle is sequencing while synthesizing. Taking the solexa sequencing method as an example, the DNA strand is randomly interrupted by physical means, and then a specific linker is added to both ends of the fragment, and an amplification primer sequence is attached to the linker.
  • DNA polymerase synthesizes the complementary strand of the fragment to be tested, and reads the base sequence by detecting the fluorescent signal carried by the newly synthesized base, thereby obtaining the sequence of the fragment to be tested (http://www.illumina.com).
  • Second-generation sequencing technology has been widely used in many fields of biological sciences, especially to study polymorphisms between different individuals in a species.
  • the traditional Call INDEL method is to compare the short reads (sequencing sequences) obtained by the sequencing individuals to the reference sequence through the comparison software, thereby obtaining the INDEL information of the sequenced individuals.
  • the available processes are: Use the SOAP software to compare the reads back to the reference sequence, and use the samtools software to process the alignment results to find the INDEL site 1 ' 2 .
  • the general process is shown in Figure 1.
  • the short sequencing fragments (or sequencing sequences) of the individual are aligned back to the reference sequence, and the base insertion of the sequencing individual relative to the reference sequence is determined by alignment (insertion of base A as shown in Figure 1). And base deletion (base T is deleted as shown in Figure 1).
  • a method of labeling a genomic single insertion deletion (INDEL) site comprising: counting the sequencing depth of a RAD single-end sequencing sequence of two individual genomes; filtering out sequencing in two individual genomes RAD single-ended sequencing sequence of depth 1; determination of heterozygous sites on RAD-tags within two individual genomes and formation of consensus sequences; non-hybrid RAD-tag sequences and heterozygous RAD-tags of two individual genomes The sequence of sequences is aligned to determine the INDEL information of the two individuals on the RAD-tag.
  • the method comprises: filtering out the RAD-tag sequence in the repeat region of the two individual genomes;
  • the alignment of the non-hybrid RAD-tag sequences of two individuals and the consensus sequence of the hybrid RAD-tag sequence comprises: non-hybrid RAD-tag sequences and heterozygous RADs in which the two individuals are in non-repetitive regions The alignment sequence of the -tag sequence is aligned.
  • the sequencing sequence has more than two copies on the genome, while the two or more copies are at different positions from the corresponding homologous chromosome There are heterozygous sites; and/or the sequencing sequence is present in multiple copies in the genome, and has a higher sequencing depth, with one copy having a heterozygous site on the corresponding homologous chromosome.
  • obtaining RAD single-ended sequencing sequences of two individual genomes by high throughput sequencing RAD single-ended sequences of the obtained two individual genomes
  • the sequencing sequence is filtered to remove unsuccessful sequencing sequences.
  • the unqualified sequencing sequence comprises: a sequencing sequence in which the number of bases whose sequencing quality is lower than a predetermined low quality threshold exceeds 50% of the number of bases of the entire sequencing sequence; and/or the sequencing result in the sequencing sequence is not a sequencing sequence in which the determined number of bases exceeds 10% of the number of bases of the entire sequencing sequence; and/or a sequencing sequence in which the exogenous sequence is present; and/or the initial few bases are not sequenced by the restriction endonuclease sequence sequence.
  • determining the heterozygous loci on the RAD-tag within the two individual genomes includes: The RAG-tags within the genome perform pairwise alignments that do not allow for gaps to determine heterozygous sites on the RAD-tag within the individual's genome.
  • the allowable mismatch number of the pairwise alignment of the unallowable voids is determined based on the length of the sequencing sequence.
  • the RAG-tag within the individual genome performs a pairwise alignment of the unallowed gaps to determine the heterozygous locus on the RAD-tag within the individual genome comprising: impermissible gaps between the sequencing sequences within the individual genome Pairwise alignment; all the sequencing sequences satisfying the alignment conditions are clustered; clustering results of only two sequencing sequences in the clustering results are selected, and the positions of the sequencing sequences are heterozygous sites.
  • the INDEL information on the RAD-tag of the two individual genomes is filtered.
  • filtering the INDEL information of the two individual genomes on the RAD-tag comprises: extracting the alignment result of the INDEL information between the two individual genomes; filtering out the INDEL number on the RAD-tag between the two individual genomes Alignment results greater than 2 and vacancy length greater than 3; for each pair of INDEL comparison results, if there are other alignment results, satisfying the alignment mismatch number is less than or equal to 3 and the vacancy number is less than or equal to 1, then The alignment result is filtered out; when the INDEL information and the enzyme cleavage site are within 3 bp of the end of the Sjf column. The result of the comparison of this type is filtered out.
  • a marker device for a genomic INDEL site comprising: a sequence depth acquisition unit for counting the sequencing depth of a RAD single-end sequencing sequence of two individual genomes; a sequence depth filtering unit, RAD single-end sequencing sequences with a sequencing depth of 1 in two individual genomes are filtered out; a consensus sequence forming unit is used to determine the heterozygous sites on the RAD-tag within the two individual genomes and form a consensus sequence; a unit for comparing the sequences of the non-hybrid RAD-tags of the two individuals and the consensus sequence of the hybrid RAD-tag to determine the INDEL information of the two individuals on the RAD-tag.
  • the apparatus further comprises: a repeating sequence filtering unit, configured to filter out RAD-tags in the repeating region inside the two individual genomes;
  • the site determining unit is a non-hybrid RAD in which the two individuals are in a non-repetitive region
  • the -tag sequence and the consensus sequence of the hybrid RAD-tag are aligned to determine the INDEL information of the two individuals on the RAD-tag.
  • the sequencing sequence has more than two copies on the genome, while the two or more copies a heterozygous site exists at a different position from the corresponding homologous chromosome; and/or the sequencing sequence is present in the genome in multiple copies, and has a higher sequencing depth, wherein one copy is heterozygous on the corresponding homologous chromosome Closing point.
  • the method further comprises: a sequencing unit for obtaining RAD single-end sequencing sequences of two individual genomes by high-throughput sequencing; a sequencing sequence filtering unit for performing RAD single-end sequencing sequences of the obtained two individual genomes Filter to remove unsuccessful sequencing sequences.
  • the unqualified sequencing sequence comprises: a sequencing sequence in which the number of bases whose sequencing quality is lower than a predetermined low quality threshold exceeds 50% of the number of bases of the entire sequencing sequence; and/or the sequencing result in the sequencing sequence is not a sequencing sequence in which the determined number of bases exceeds 10% of the number of bases of the entire sequencing sequence; and/or a sequencing sequence in which the exogenous sequence is present; and/or the initial few bases are not sequenced by the restriction endonuclease sequence sequence.
  • the consensus sequence forming unit performs a pairwise alignment of the RAG-tags within the individual genome with no allowable gaps to determine the heterozygous sites on the RAD-tag within the individual genome and form a consensus sequence.
  • the number of mismatches allowed for pairwise alignment of the gaps is determined based on the length of the sequencing sequence.
  • the apparatus further comprises: an INDEL information filtering unit for filtering INDEL information of the two individual genomes on the RAD-tag.
  • the INDEL information filtering unit extracts the alignment result of the INDEL information between the two individual genomes; filters out the alignment result of the INDEL number on the RAD-tag between the two individual genomes greater than 2 and the vacancy length greater than 3; The alignment result of each pair of INDEL, if there are other alignment results, if the mismatch number of the matched pair is less than or equal to 3 and the gap number is less than or equal to 1, the alignment result is filtered out; when the INDEL information and the enzyme are cut The distance between the site and the end of the sequence is within 3 bp. Then the comparison result of this type is filtered out.
  • the method and device for labeling the genomic INDEL locus in the embodiment of the present invention realizes a bioinformatics analysis method for accurately searching for INDEL locus markers in a certain species population by processing RAD sequencing data in the absence of a reference sequence.
  • FIG. 1 is a schematic view showing a method of determining a SNP site in the prior art
  • 2A-2D are schematic diagrams showing various steps of a RAD sequencing technique
  • Figure 3 is a schematic diagram showing an example of RAD single-ended sequencing of a genome
  • Figure 4 is a flow chart showing one embodiment of the genomic INDEL site labeling method of the present invention.
  • Figure 5 is a flow chart showing another embodiment of the genomic INDEL site labeling method of the present invention.
  • Figure 6 is a schematic diagram showing the statistical information of the depth information of the sequencing sequence
  • Figure 7 is a schematic diagram showing the depth information storage of a sequencing sequence
  • Figure 8 is a diagram showing an example of determining hybrid site information inside an individual based on the first drawer principle
  • Figure 9 is a schematic view showing an example of a heterozygous site located in a repeating region
  • Figure 10 shows a RAD-tag sequencing depth profile
  • Figure 11 is a view showing the structure of one embodiment of the genomic INDEL site labeling device of the present invention.
  • Fig. 12 is a view showing the configuration of another embodiment of the genomic INDEL site labeling device of the present invention. detailed description
  • the inventors of the present application developed a bioinformatics analysis method to process RAD (Restriction-site Associated DNA) data to find RAD fragments.
  • RAD Restriction-site Associated DNA
  • the INDEL locus information simplifies the complexity of the genome and reduces the cost of sequencing 3 .
  • FIGS. 2A–2D show schematic diagrams of the various steps of RAD sequencing technology.
  • the experimental procedure for the development of RAD molecular markers includes: P-type endonuclease digestion of genomic DNA, plus PI linker, P1 linker containing primer sequences required for amplification, Illumina sequencing primer binding sites Sequences and short tag sequences that distinguish different samples (Figure 2A); samples with different P1 linkers are mixed together, physically broken into 300-700 bp sequences ( Figure 2B); P2 linker added ( Figure 2C); PCR Amplify and enrich RAD tags (Fig. 2D) to construct a library of high-throughput sequencing.
  • a hash table (also called a hash table) is a data structure that is accessed directly based on a key value. That is, it accesses the record by mapping the key value to a location in the table to speed up the lookup.
  • This mapping function is called a hash function, and the array in which the records are stored is called a hash table. Indexing data using a hash table basically grows linearly as the amount of data increases, and a string of "ATCGN" has a very low probability of conflicting key values. This has good performance when dealing with massive sequencing data.
  • n drawers there are at least one drawer with 2 or more objects. From this principle we can deduce that if n-1 objects are placed in n drawers, there is at least one object in the drawer.
  • the bioinformatics analysis method of the INDEL site information on the inter-RAD segment is designed to overcome some of the technical bottlenecks of the traditionally obtained INDEL method.
  • FIG. 3 shows the palindrome sequence of the "G A AATTC” on the DNA molecule identified by the restriction endonuclease Ecorl, and the DNA molecule is cleaved between G and VIII, and the enzyme is cleaved.
  • the DNA molecule is physically broken into short sequence fragments, and a ligation end is added to the end of the restriction enzyme, and the DNA fragment is single-end sequenced.
  • the sequencing read length is generally 50 nt or 100 nt.
  • Enrichment and sequencing of specific regions of the genome will be done by RAD sequencing, which will reduce the complexity of the genome and the cost of sequencing.
  • Figure 4 is a flow chart showing one embodiment of the genomic INDEL site labeling method of the present invention.
  • step 402 the sequencing depth of the RAD single-end sequencing sequences of the two individual genomes is counted.
  • the sequencing depth of the RAD single-end sequencing sequence for each individual genome was counted.
  • Step 404 filtering out the RAD single-end sequencing sequence with a sequencing depth of 1 in the two individual genomes.
  • short sequences of depth 1 are caused by sequencing errors.
  • short sequence information of depth 1 is filtered out, reducing false detection results due to sequencing errors.
  • step 406 the heterozygous sites on the RAD-tag within the two individual genomes are determined and a consensus sequence is formed.
  • the hybrid RAD-tag "GAATTCACCC and "GAATTCACIC” in the individual genome will be represented as the consensus sequence "GAATTCACSC". It is used to indicate that the individual has a heterozygous site at that position.
  • Step 408 Align the non-hybrid RAD-tag sequences of the two individual genomes with the consensus sequence of the hybrid RAD-tag sequence to determine the INDEL information of the two individuals on the RAD-tag. This alignment refers to the alignment between the sequences of two individual genomes.
  • the bioinformatics analysis method for accurately detecting the INDEL locus in a certain species population by using the RAD sequencing data in the absence of the reference sequence is realized, and some technical bottlenecks of the conventional INDEL method are overcome.
  • Figure 5 shows another embodiment of the genomic INDEL site labeling method of the present invention. Flow chart.
  • step 502 RAD single-end sequencing sequences of two individual genomes are obtained by high-throughput sequencing, wherein the high-throughput sequencing technology can be Illumina GA sequencing technology, or other existing high-throughput sequencing. technology.
  • Step 504 filtering the obtained RAD single-end sequencing sequences of the two individual genomes to remove the unqualified sequencing sequences.
  • the sequencing sequence is filtered to remove the unqualified sequence.
  • the unqualified sequence includes, for example, that the number of bases whose sequencing quality is below a certain threshold (e.g., the single base sequencing quality is less than 20) exceeds 50% of the number of bases of the entire sequence, and is considered to be a non-conforming sequence.
  • the low quality threshold is determined by the specific sequencing technology and the sequencing environment; the number of bases with undetermined sequencing results in the sequence (such as N in Illumina GA sequencing results) exceeds 10% of the number of bases in the entire sequence and is considered unqualified.
  • sequence In addition to the sample linker sequence, it is aligned with other experimentally introduced exogenous sequences, such as various linker sequences. If the exogenous sequence is present in the sequence, it is considered to be a non-conforming sequence; in the sequence, if the first few bases are not the end-cut sequence, it is filtered out (such as the restriction endonuclease Ecorl, if the short sequence starts at the beginning) AATTC" filters out the entire short sequence).
  • exogenous sequence is present in the sequence, it is considered to be a non-conforming sequence; in the sequence, if the first few bases are not the end-cut sequence, it is filtered out (such as the restriction endonuclease Ecorl, if the short sequence starts at the beginning) AATTC" filters out the entire short sequence).
  • Step 506 counting the sequencing depth of the RAD single-end sequencing sequences of the two individual genomes. The same short sequences in each individual genome were counted separately.
  • the short sequence information filtered by the individual is used as a hash key, and the short sequence is counted by the hash value.
  • the hash value can be implemented in any programming language, such as a C++ hash table) so that you can get the sequencing depth information for each short sequence in an individual. The specific process is shown in Figure 6.
  • the information of the heap is saved as shown in Fig. 7.
  • the first column indicates the RAD sequence information; the second column indicates the number of times the sequence is sequenced, that is, the depth information; and the third column indicates the ID of the sequence information.
  • Step 508 filtering out the RAD single-end sequencing sequence with a sequencing depth of 1 in the two individual genomes.
  • step 510 the heterozygous sites on the RAD-tag within the two individual genomes are determined and a consensus sequence is formed.
  • the RAD-tag inside the individual genome performs a pairwise alignment of the vacancies, looking for information on the heterozygous loci within the individual genome. Generally, if the sequencing length is less than 50 nt, the allowed error The match number is 1.
  • the short sequence in the individual is segmented, and a hash table is built according to the divided substring to perform indexing. If one mismatch is allowed, the short sequence of one individual is averaged into two substrings, so that a short sequence can have a mismatch ratio with another short sequence in the individual. According to the drawer principle, the mismatch is either On the left, or on the right, there must be one side where there is no mismatch ratio. That is to say, if m mismatches are allowed, they are divided into m+1 substrings, then at least one substring is not mismatched and can be completely compared. In this case, you can use the split substring as a seed to create a hash table.
  • the average divided substring is used as the hash key, and the entire string is used as the hash value to create a hash table to index the string.
  • most of the strings close to the string can be quickly found through the hash table, narrowed down by the hash table, and then aligned one by one to find the heterozygous site information inside the individual.
  • the specific process is shown in Figure 8. (The string alignment algorithm can be implemented in any programming language, such as C++).
  • step 512 the RAD-tag sequence within the repeating region of the two individual genomes is filtered out.
  • the RAD-tags that can be compared between the two internal groups are clustered, and the clustering results of only one RAD-tag and two RAD-tags in the clustering result are selected. Only one RAD-tag clustering result indicates that there is no heterozygous locus in the sequencing region, and only two RAD-tag clustering results indicate that there is a heterozygous locus in the sequencing region, and in general, it will be in the genome. Non-repeating area.
  • (b) shows that sequence 1 has multiple copies on the genome, one of which has a heterozygous site on the corresponding homologous chromosome, and the alignment result of (b) appears when the alignment is performed.
  • the more complicated cases caused by other repetitive sequences are based on these two cases, and the RAD-tag of the repeated area is filtered out during the processing.
  • Step 514 Align the non-hybrid RAD-tag sequences of the two individuals in the non-repetitive region with the consensus sequence of the hybrid RAD-tag sequence to determine the INDEL information of the two individuals on the RAD-tag.
  • a lookup of the INDEL locus on the RAD-tag between two sequencing individuals Through the processing in the previous step, the heterozygous RAD-tag in the body will be shown as a consistent sequence, and the information of the repeated sequence will be filtered out.
  • the sequence of the non-hybrid RAD-tag in which the two individuals are in the non-repetitive region is aligned with the consensus sequence of the hybrid RAD-tag to find the INDEL information of the two individuals on the RAD-tag, using the comparison software
  • the parameters of the comparison software are generally default.
  • Step 516 filtering the INDEL information of the two individual genomes on the RAD-tag
  • blat is used for comparison, and the results are filtered to find high-reliability INDEL locus information.
  • the 50bp sequencing length of the RAD-tag process is as follows:
  • the number of mismatches to be matched is less than or equal to 3 and the number of vacancies is less than or equal to 1. Analyze whether there are other alignment results, and if so, prove that the RAD-tag is repeated to multiple locations and the alignment is filtered out. Because such alignment results are likely due to the short sequence of sequencing and the interference of repeated sequence information, the alignment bias caused.
  • the accuracy and reliability of the INDEL locus information detection are improved by filtering the sequencing sequence, excluding the RAD-tag sequence located in the repetitive region, and filtering the INDEL information on the RAD-tag.
  • Fig. 11 is a view showing the configuration of an embodiment of the genomic INDEL site labeling device of the present invention.
  • the apparatus includes: a sequence depth obtaining unit 113 for counting the sequencing depth of RAD single-end sequencing sequences of two individual genomes; and a sequence depth filtering unit 115 for filtering out sequencing sequences with a sequencing depth of 1
  • the consensus sequence forming unit 117 determines the heterozygous sites on the RAD-tags within the two individual genomes and forms a consensus sequence. For example, the consensus sequence forming unit 117 performs a pairwise alignment of the RAG-tags within the individual genome that does not allow gaps to determine the heterozygous sites on the RAD-tag within the individual genome and form a consensus sequence.
  • a site determining unit 119 is configured to compare the sequence of the non-hybrid RAD-tag of the two individuals with the consensus sequence of the hybrid RAD-tag to determine the INDEL information of the two individuals on the RAD-tag.
  • Fig. 12 is a view showing the configuration of another embodiment of the genomic INDEL site marker device of the present invention.
  • the INDEL site marking device includes a sequence depth obtaining unit 113, a sequence depth filtering unit 115, a consistency sequence forming unit 117, a repetition sequence filtering unit 128, and a site determining unit 129.
  • the sequence depth obtaining unit 113, the sequence depth filtering unit 115, and the consistency sequence forming unit 117 can refer to the description of the foregoing embodiment, and will not be described in detail herein for the sake of brevity.
  • the repeat sequence filtering unit 128 is configured to filter out RAD-tags within the repeating regions of the two individual genomes.
  • the sequencing sequence has more than two copies on the genome, while the two or more copies are heterozygous at different positions with the corresponding homologous chromosomes The binding site; and/or the sequencing sequence is present in the genome in multiple copies with a high depth of sequencing, one of which is stored on the corresponding homologous chromosome At the heterozygous site.
  • the site determining unit 129 aligns the non-hybrid RAD-tag sequences of the two individuals in the non-repetitive region with the consensus sequence of the hybrid RAD-tag.
  • the INDEL site labeling device further comprises: a sequencing unit 121 for obtaining RAD single-ended sequencing sequences of two individual genomes; a sequencing sequence filtering unit 122 for RAD singles of the obtained two individual genomes
  • the end sequencing sequence is filtered to remove unqualified sequencing sequences.
  • the unqualified sequencing sequence includes, for example, a sequencing sequence in which the number of bases whose sequencing quality is lower than a predetermined low quality threshold exceeds 50% of the number of bases of the entire sequencing sequence; and/or the base in which the sequencing result is indeterminate in the sequencing sequence. a sequencing sequence having a number exceeding 10% of the number of bases of the entire sequencing sequence; and/or a sequencing sequence in which the exogenous sequence is present; and/or a starting sequence of several bases that is not a restriction endonuclease sequence.
  • the apparatus further includes an INDEL information filtering unit 130 for filtering INDEL information on the RAD-tag of the two individual genomes.
  • the INDEL information filtering unit 130 extracts the alignment result with the INDEL information between the two individual genomes; filters out the alignment result of the INDEL number on the RAD-tag between the two individual genomes greater than 2 and the gap length greater than 3; The alignment result of a pair of INDELs, if there are other alignment results, if the mismatch number of the matched pair is less than or equal to 3 and the gap number is less than or equal to 1, the alignment result is filtered out; when the INDEL information and the enzyme cleavage position The distance at the end of the SJ ⁇ column is within 3 bp. Then the comparison result of this type is filtered out.

Landscapes

  • Chemical & Material Sciences (AREA)
  • Organic Chemistry (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Zoology (AREA)
  • Proteomics, Peptides & Aminoacids (AREA)
  • Health & Medical Sciences (AREA)
  • Engineering & Computer Science (AREA)
  • Wood Science & Technology (AREA)
  • Analytical Chemistry (AREA)
  • Microbiology (AREA)
  • Physics & Mathematics (AREA)
  • Molecular Biology (AREA)
  • Immunology (AREA)
  • Biotechnology (AREA)
  • Biophysics (AREA)
  • Biochemistry (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • General Engineering & Computer Science (AREA)
  • General Health & Medical Sciences (AREA)
  • Genetics & Genomics (AREA)
  • Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)

Abstract

本发明公开了一种基因组INDEL位点的标记方法，其包括：统计两个个体基因组的RAD单端测序序列的测序深度，过滤掉两个个体基因组中测序深度为1的RAD单端序列，确定两个个体基因组内部RAD-标记上的杂合位点并形成一致性序列，将两个个体基因组的非杂合RAD-标记序列及杂合RAD-标记序列的一致性序列进行比对，以确定两个个体在RAD-标记上的INDEL信息。本发明还公开了一种基因组INDEL位点的标记装置。

Description

基因组 I NDEL位点标记方法和装置 技术领域
本发明属于基因工程技术领域， 尤其涉及一种基因组 INDEL ( Insertion -Deletion, ***缺失)位点标记方法和装置。 背景技术
基因组 INDEL位点标记， 指的是两种亲本中在全基因组中的差异， 相对另一个亲本而言， 其中一个亲本的基因组中有一定数量的核苷酸插 入或缺失。 INDEL 位点信息的获得可以有许多重要的应用， 如构建遗传 图谱， 基因分型， 分子标记育种， 疾病检测等。
如今， 第二代 DNA测序技术是一种高通量低成本的测序技术， 基本 原理是边合成边测序。 以 solexa测序方法为例， 先用物理方法将 DNA链 随机打断， 然后在片段两端加上特定接头， 接头上有扩增引物序列。 测 序时， DNA 聚合酶合成待测片段的互补链， 通过检测新合成碱基所携带 的 荧光信号读取碱基序列 ， 从而获得待测 片 段的序列 ( http://www.illumina.com )。
第二代测序技术已经广泛应用于生物科学的许多领域， 特别是研究 一个物种不同个体之间的多态性。 传统 Call INDEL的方法是将测序个体 得到的短 reads (测序序列）通过比对软件比对回参考序列， 从而得到测 序个体的 INDEL信息。 可用的流程有： 使用 SOAP软件将 reads比对回 参考序列， 使用 samtools软件处理比对结果寻找 INDEL位点 1 '2。 大体过 程如图 1所示， 将个体短的测序片段（或测序序列） 比对回参考序列， 通 过比对确定测序个体相对于参考序列的碱基***（如图 1 所示***碱基 A )和碱基缺失（如图 1所示缺失碱基 T )。
目前， 有参考序列的物种都可以很方便的进行 INDEL标记的开发， 但是对于那些非模式生物基本上是没有参考序列的。 在没有参考序列的 情况下， 传统获得 INDEL的方法存在着技术上的瓶颈。 发明内容 本发明的发明人发现上述现有技术中存在问题， 并因此针对所述 问题中的至少一个问题提出了一种新的技术方案。
本发明的一个目的是提供一种用于基因组 INDEL位点标记的技术 方案。
根据本发明的第一方面， 提供了一种基因组单***缺失（INDEL ) 位点的标记方法， 包括： 统计两个个体基因组的 RAD单端测序序列的测 序深度； 过滤掉两个个体基因组中测序深度为 1的 RAD单端测序序列； 确定两个个体基因组内部 RAD-tag 上的杂合位点并形成一致性序列； 将 两个个体基因组的非杂合 RAD-tag序列及杂合 RAD-tag序列的一致性序 列进行比对， 以确定两个个体在 RAD-tag上的 INDEL信息。
可选地， 在将两个个体的非杂合 RAD-tag 的序列及杂合 RAD-tag 的 一致性序列进行比对之前包括： 过滤掉两个个体基因组内部处于重复区 域的 RAD-tag序列；所述将两个个体的非杂合 RAD-tag序列及杂合 RAD- tag序列的一致性序列进行比对包括： 将两个个体处于非重复区域的非杂 合 RAD-tag序列及杂合 RAD-tag序列的一致性序列进行比对。
可选地， 满足如下条件作为处于基因组序列的重复区域中的 RAD-tag 序列： 测序序列在基因组上存在两个以上拷贝， 同时所述两个以上拷贝 与对应的同源染色体在不同的位置上存在杂合位点； 和 /或测序序列在基 因组在存在多个拷贝， 且具有较高的测序深度， 其中一个拷贝与对应的 同源染色体上存在杂合位点。
可选地， 在统计两个个体基因组的 RAD单端测序序列的测序深度之 前包括： 通过高通量测序获得两个个体基因组的 RAD单端测序序列； 对 获得的两个个体基因组的 RAD单端测序序列进行过滤以去除不合格的测 序序列。
可选地， 不合格的测序序列包括： 测序质量低于预定的低质量阈值的 碱基个数超过整条测序序列碱基个数的 50%的测序序列； 和 /或测序序列 中测序结果不确定的碱基个数超过整条测序序列碱基个数的 10%的测序 序列； 和 /或存在外源序列的测序序列； 和 /或起始的几个碱基不是酶切末 端序列的测序序列。
可选地， 确定两个个体基因组内部 RAD-tag上的杂合位点包括： 个体 基因组内部的 RAG-tag进行不容许空隙的两两比对以确定个体基因组内 部 RAD-tag上的杂合位点。
可选地， 根据测序序列的长度确定所述不容许空隙的两两比对的容许 的错配数。
可选地， 个体基因组内部的 RAG-tag进行不容许空隙的两两比对以确 定个体基因组内部 RAD-tag 上的杂合位点包括： 在个体基因组内部的测 序序列之间进行不容许空隙的两两比对； 将所有满足比对条件的测序序 列进行聚类； 挑选出聚类结果中只有两种测序序列的聚类结果， 该测序 序列的位置即存在杂合位点。
可选地， 对两个个体基因组在 RAD-tag上的 INDEL信息进行过滤。 可选地， 对两个个体基因组在 RAD-tag上的 INDEL信息进行过滤包 括： 提取两个个体基因组之间具有 INDEL信息的比对结果； 过滤掉两个 个体基因组之间 RAD-tag上 INDEL数大于 2以及空位长度大于 3的比对 结果； 对每一对的 INDEL 的比对结果， 如果存在其他的比对结果， 满足 比对的错配数小于等于 3 并且空位数小于等于 1 , 则将该比对结果过滤 掉； 当 INDEL信息与酶切位点以 Sjf列末端的距离在 3bp以内。 则将该 类型的比对结果过滤掉。
根据本发明的另一方面， 提供一种基因组 INDEL 位点的标记装置， 包括： 序列深度获取单元， 用于统计两个个体基因组的 RAD单端测序序 列的测序深度； 序列深度过滤单元， 用于过滤掉两个个体基因组中测序 深度为 1的 RAD单端测序序列； 一致性序列形成单元， 用于确定两个个 体基因组内部 RAD-tag 上的杂合位点并形成一致性序列； 位点确定单 元， 用于将两个个体的非杂合 RAD-tag的序列及杂合 RAD-tag的一致性 序列进行比对以确定两个个体在 RAD-tag上的 INDEL信息。
可选地， 该装置还包括： 重复序列过滤单元， 用于过滤掉两个个体基 因组内部处于重复区域的 RAD-tag;所述位点确定单元将两个个体处于非 重复区域的非杂合 RAD-tag序列及杂合 RAD-tag的一致性序列进行比对 以确定两个个体在 RAD-tag上的 INDEL信息。
可选地， 满足如下条件作为处于基因组序列的重复区域中的 RAD-tag 序列： 测序序列在基因组上存在两个以上拷贝， 同时所述两个以上拷贝 与对应的同源染色体在不同的位置上存在杂合位点； 和 /或测序序列在基 因组在存在多个拷贝， 且具有较高的测序深度， 其中一个拷贝与对应的 同源染色体上存在杂合位点。
可选地， 还包括： 测序单元， 用于通过高通量测序获得两个个体基因 组的 RAD单端测序序列； 测序序列过滤单元， 用于对获得的两个个体基 因组的 RAD单端测序序列进行过滤以去除不合格的测序序列。
可选地， 不合格的测序序列包括： 测序质量低于预定的低质量阈值的 碱基个数超过整条测序序列碱基个数的 50%的测序序列； 和 /或测序序列 中测序结果不确定的碱基个数超过整条测序序列碱基个数的 10%的测序 序列； 和 /或存在外源序列的测序序列； 和 /或起始的几个碱基不是酶切末 端序列的测序序列。
可选地， 一致性序列形成单元对个体基因组内部的 RAG-tag进行不容 许空隙的两两比对以确定个体基因组内部 RAD-tag 上的杂合位点并形成 一致性序列。 不容许空隙的两两比对的容许的错配数根据测序序列的长 度确定。
可选地， 该装置还包括： INDEL信息过滤单元， 用于对两个个体基因 组在 RAD-tag上的 INDEL信息进行过滤。
可选地， INDEL 信息过滤单元提取两个个体基因组之间具有 INDEL 信息的比对结果； 过滤掉两个个体基因组之间 RAD-tag上 INDEL数大于 2以及空位长度大于 3的比对结果； 对每一对的 INDEL的比对结果， 如 果存在其他的比对结果， 满足比对的错配数小于等于 3并且空位数小于等 于 1 , 则将该比对结果过滤掉； 当 INDEL信息与酶切位点以及序列末端 的距离在 3bp以内。 则将该类型的比对结果过滤掉。
本发明实施例中的基因组 INDEL 位点的标记方法和装置， 实现了缺 少参考序列的情况下， 处理 RAD 测序数据准确寻找某个物种群体中 INDEL位点标记的生物信息学分析方法。
通过以下参照附图对本发明的示例性实施例的详细描述， 本发明 的其它特征及其优点将会变得清楚。 附图说明 构成说明书的一部分的附图描述了本发明的实施例， 并且连同说 明书一起用于解释本发明的原理。
参照附图， 根据下面的详细描述， 可以更加清楚地理解本发明， 其中：
图 1示出现有技术中的 SNP位点的确定方法的示意图；
图 2A~2D示出 RAD测序技术的各个步骤的示意图；
图 3示出基因组的 RAD单端测序的一个例子的示意图；
图 4示出本发明的基因组 INDEL位点标记方法的一个实施例的流 程图；
图 5示出本发明的基因组 INDEL位点标记方法的另一个实施例的 流程图；
图 6示出测序序列的深度信息统计示意图；
图 7示出测序序列的深度信息存储示意图；
图 8示出基于第一抽屉原理确定个体内部的杂合位点信息的例子的 示意图；
图 9示出位于重复区域的杂合位点的例子的示意图；
图 10示出 RAD-tag测序深度分布图；
图 11 示出本发明的基因组 INDEL位点标记装置的一个实施例的 结构图；
图 12示出本发明的基因组 INDEL位点标记装置的另一个实施例 的结构图。 具体实施方式
现在将参照附图来详细描述本发明的各种示例性实施例。 应注意 到： 除非另外具体说明， 否则在这些实施例中阐述的部件和步骤的相对 布置、 数字表达式和数值不限制本发明的范围。
同时， 应当明白， 为了便于描述， 附图中所示出的各个部分的尺寸 并不是按照实际的比例关系绘制的。
以下对至少一个示例性实施例的描述实际上仅仅是说明性的， 决不 作为对本发明及其应用或使用的任何限制。 对于相关领域普通技术人员已知的技术、 方法和设备可能不作详细 讨论， 但在适当情况下， 所述技术、 方法和设备应当被视为授权说明书 的一部分。
在这里示出和讨论的所有示例中， 任何具体值应被解释为仅仅是示 例性的， 而不是作为限制。 因此， 示例性实施例的其它示例可以具有不 同的值。
应注意到： 相似的标号和字母在下面的附图中表示类似项， 因此， 一旦某一项在一个附图中被定义， 则在随后的附图中不需要对其进行进 一步讨论。
为了突破非模式生物缺少参考序列的瓶颈， 本申请发明人开发了一 套生物信息学分析方法， 处理 RAD ( Restriction-site Associated DNA, 限制性内切位点相关 DNA )数据， 寻找 RAD 片段上的 INDEL 位点信 息， 简化了基因组的复杂度， 同时也减少了测序成本 3。
RAD测序技术采用了新的建库方式， 图 2A~2D示出 RAD测序技术 的各个步骤的示意图。 如图 2A~2D所示， RAD分子标记开发的实验流程 包括： P艮制性内切酶消化基因组 DNA， 加上 PI接头， P1接头包含扩增 所需的引物序列、 Illumina测序引物结合位点序列以及区分不同样品的短 标签序列 （图 2A ); 带有不同 P1接头的样品混合在一起， 物理方法打断 成 300-700bp 的序列 （图 2B ); 加上 P2接头（图 2C ); PCR扩增富集 RAD tags (图 2D )， 从而构建上机文库进行高通量测序。
散列表（Hash table, 也叫哈希表）， 是根据关键码值 (Key value)而直 接进行访问的数据结构。 也就是说， 它通过把关键码值映射到表中一个 位置来访问记录， 以加快查找的速度。 这个映射函数叫做散列函数， 存 放记录的数组叫做散列表。 使用哈希表对数据进行索引基本是随着数据 量的上升线性增长， 而且由 "ATCGN" 构成的字符串， 键值出现冲突的 可能性非常低。 这样在处理海量测序数据的时候有着很好的性能。
第一抽屉原理， 把多于 n 个的物体放到 n 个抽屉里， 则至少有一个 抽屉里有 2个或 2个以上的物体。 由这个原理我们可以推导出， 如果把 n-1个的物体放到 n个抽屉中， 则至少有一个抽屉里是没有物体。
本发明的目的在于提供一种直接处理 RAD测序数据， 寻找两个个体 间 RAD片段上的 INDEL位点信息的生物信息学分析方法， 旨在克服传 统获得 INDEL方法的一些技术瓶颈。
基于先期完整的实验流程， 将得到两个个体基因组的 RAD单端测序 数据。 如图 3所示， 该图显示了用限制性内切酶 Ecor l, 识别 DNA分子 上 "GAAATTC" 的回文序列， 并在 G与八之间将 DNA分子切断， 将酶 切后的 DNA分子用物理方法打断成短的序列片段， 并在其中酶切的一端 加上接头并对 DNA片段进行单末端测序， 测序读长一般为 50nt, 也可以 为 100nt。
通过 RAD测序方式将会对基因组的特定区域进行富集测序， 这样做 会降低基因组的复杂度和测序的成本。
图 4示出本发明的基因组 INDEL位点标记方法的一个实施例的流 程图。
如图 4所示， 步骤 402， 统计两个个体基因组的 RAD单端测序序列 的测序深度。 统计每个个体基因组的 RAD单端测序序列的测序深度。
步骤 404, 过滤掉两个个体基因组中测序深度为 1的 RAD单端测序 序列。 一般来讲， 深度为 1的短序列是由测序错误导致的， 在这一步过滤 掉深度为 1的短序列信息， 减少由于测序错误引起的错误的检测结果。
步骤 406， 确定两个个体基因组内部 RAD-tag 上的杂合位点并形成 一致性序列。 对两个测序个体内部 RAD-tag ( Restriction -site Associated DNA tag, P艮制性内切位点相关 DNA标记）上杂合位点信息进行整合， 生成一致性序列 。 比如， 个体基因 组 内 的杂合 RAD-tag " GAATTCACCC 和 " GAATTCACIC " 将被表示成一致性序列 "GAATTCACSC"。 用于表示个体在该位置存在一个杂合位点。
步骤 408, 将两个个体基因组的非杂合 RAD-tag序列及杂合 RAD- tag 序列的一致性序列进行比对， 以确定两个个体在 RAD-tag 上的 INDEL信息。 该比对指两个个体基因组的序列之间的比对。
上述实施例中， 实现了缺少参考序列的情况下， 处理 RAD测序数据 准确寻找某个物种群体中 INDEL位点标记的生物信息学分析方法， 克服 传统获得 INDEL方法的一些技术瓶颈。
图 5示出本发明的基因组 INDEL位点标记方法的另一个实施例的 流程图。
如图 5所示， 步骤 502, 通过高通量测序获得两个个体基因组的 RAD 单端测序序列， 其中高通量测序技术可以为 Illumina GA测序技术， 也可 以为现有的其他高通量测序技术。
步骤 504, 对获得的两个个体基因组的 RAD单端测序序列进行过滤 以去除不合格的测序序列。 接收到高通量测序序列后， 对测序序列进行 过滤， 去除不合格的序列。 不合格序列例如包括： 测序质量低于某一阈 值（如单碱基测序质量低于 20 ) 的碱基个数超过整条序列碱基个数的 50%则认为是不合格序列。 低质量阈值由具体测序技术及测序环境而定； 序列中测序结果不确定的碱基（如 Illumina GA测序结果中的 N )个数超 过整条序列碱基个数的 10%则认为是不合格序列； 除样本接头序列外， 与其它实验引入的外源序列比对， 如各种接头序列。 若序列中存在外源 序列则认为是不合格序列； 在序列中， 若起始的几个碱基不是酶切末端 序列则过滤掉（如限制性内切酶 Ecor l , 短序列开头若不是 "AATTC" 则过滤掉整个短序列）。
步骤 506, 统计两个个体基因组的 RAD单端测序序列的测序深度。 分别对每个个体基因组中相同的短序列进行统计计数。
在本发明的一个实施例中， 将个体过滤后的短序列信息作为哈希的 键， 用哈希的值对短序列进行计数。 （可以用任何一种编程语言实现， 如 C++的哈希表)这样就可以得到一个个体中每一种短序列的测序深度信息。 具体过程如图 6所示。
堆的信息以图 7的方式保存， 第一列表示的是 RAD序列信息； 第二 列表示的是该序列被测序的次数， 即深度信息； 第三列是该序列信息的 ID。
步骤 508, 过滤掉两个个体基因组中测序深度为 1的 RAD单端测序 序列。
步骤 510， 确定两个个体基因组内部 RAD-tag 上的杂合位点并形成 一致性序列。
个体基因组内部的 RAD-tag 进行不开空位的两两比对， 寻找个体基 因组内部的杂合位点信息。 一般测序长度小于 50nt 的情况下， 容许的错 配数为 1。
对个体内的短序列进行分割 , 并根据分割后的子字符串建立哈希表 进行索引。 如果容许 1个错配， 就将其中一个个体的短序列平均切成两个 子串， 这样某个短序列与个体内另一个短序列能够存在一个错配比上的 话， 根据抽屉原理， 错配要么在左边， 要么在右边， 这样肯定有一边是 不存在错配比上的地方。 也就是说， 如果容许 m个错配， 就分割成 m+1 个子串， 那么至少有一个子串是不存在错配能够完全比上。 这样的话， 可以将分割后的子字符串作为种子， 建立哈希表。 比如， 容许一个错配 的话， 就用平均分割后的子串作为哈希的键， 整个字符串作为哈希的 值， 建立一个哈希表， 实现对字符串的索引。 这样在处理字符串比对的 时候可以迅速通过哈希表找到大部分与该字符串相近的字符串 , 通过哈 希表缩小范围后再逐一进行比对， 找到个体内部的杂合位点信息。 具体 过程如图 8 所示 (该字符串比对算法可以通过任何一种编程语言实现， 比 如 C++)。
步骤 512， 过滤掉两个个体基因组内部处于重复区域的 RAD-tag序 列。
将个体内部两两之间能够比对上的 RAD-tag 进行聚类， 挑选出聚类 结果中只有一条 RAD-tag 和两条 RAD-tag 的聚类结果。 其中只有一条 RAD-tag 的聚类结果表明在测序区域不存在杂合位点， 只有两条 RAD- tag 的聚类结果表明在测序区域存在存在杂合位点， 且在一般情况下会处 于基因组的非重复区域。
通常， 存在如图 9所示的聚类结果及其衍生结果， 就认为这些聚类结 果中的所有 RAD-tag是处于基因组 DNA序列的重复区域， 通常会把这些 RAD-tag过滤掉。
( a )显示的是， 序列 2在基因组上存在两个拷贝， 同时这两个拷贝 与对应的同源染色体， 在不同的位置上存在杂合位点。 这样， 在比对聚 类的时候就会出现 ) 的比对结果。
( b )显示的是， 序列 1 在基因组上存在多个拷贝， 其中一个拷贝上 与对应的同源染色体上存在杂合位点， 比对的时候就会出现（b ) 的比对 结果。 其他重复序列导致的更复杂的情况， 都是以这两种情况作为基础 的， 在处理过程中都会把重复区域的 RAD-tag过滤掉。
在过滤处于重复区域的 RAD-tag 之后， 需要对具有杂合位点的 RAD-tag进行信息的整合， 生成一致性序列。
步骤 514， 将两个个体处于非重复区域的非杂合 RAD-tag序列及杂 合 RAD-tag序列的一致性序列进行比对以确定两个个体在 RAD-tag上的 INDEL信息。
两个测序个体之间 RAD-tag上 INDEL位点的查找。 通过上一步的处 理， 个体内的杂合 RAD-tag将会^ ^示成一致性序列， 重复序列的信息 将会被过滤掉。
将两个个体处于非重复区域的非杂合 RAD-tag 的序列和杂合 RAD- tag的一致性序列进行比对， 以寻找两个个体在 RAD-tag上的 INDEL信 息， 使用的比对软件可以是任何一款序列比对软件， 如 blast、 blat. 比对 软件的参数一般情况下都是默认。
步骤 516， 对两个个体基因组在 RAD-tag上的 INDEL 信息进行过 滤
例如， 采用 blat 进行比对， 对比对结果进行过滤处理以寻找高可信 度的 INDEL位点信息。
通常情况下， 50bp测序长度的 RAD-tag处理过程如下所示：
( 1 )首先需要把两个个体之间具有 INDEL 信息的比对结果的信息 提取出来。
( 2 )过滤掉两个个体之间 RAD-tag上 INDEL数大于 2以及空位长 度大于 3的比对结果。
( 3 )对每一对的 INDEL 的比对结果， 需要满足比对的错配数小于 等于 3并且空位数小于等于 1。 分析是否存在其他的比对结果， 如果存在 的话， 就证明该 RAD-tag 重复比对到多个地方， 就将该比对结果过滤 掉。 因为这样的比对结果很可能是由于测序的序列过短， 以及重复序列 信息的干扰， 造成的比对偏差。
( 4 )接下去查看过滤之后的 INDEL 距离酶切位点的距离以及序列 末端的距离。 如果 INDEL与酶切位点以 SJ^列末端的距离在 3bp以内。 则将该类型的比对结果过滤掉。
通过上述两个个体 RAD-tag数据的内部比对， 聚类， 重复区域的筛 选， 以及比对结果的过滤。 最终会得到两个个体之间具有足够深度信息 支持的 RAD-tag INDEL标记集合。
上述实施例中， 通过对测序序列进行过滤、 排除位于重复区域的 RAD-tag序列、 对 RAD-tag上的 INDEL信息进行过滤等步骤， 提高了 INDEL位点信息检测的准确度和可靠性。
本领域技术任意应当理解， 上述实施例中的步骤 502、 504、 512、 516等步骤可以是可选步骤， 或者在不同的实施例中包括其中一个或者多 个步骤。 下面介绍根据本发明基因组 INDEL位点标记方法的一个应用例。 实施例数据：
羽扇豆自交群体两个亲本的 RAD-tag测序数据。
实施例具体操作¾½：
1 )将两个亲本 RAD-tag的测序数据， 根据测序质量值， N的含量， 以及是否含有酶切末端序列进行过滤， 去除不合格的序列， 最后得到的 有效数据统计如表 1所示：
Figure imgf000013_0001
表 1、 羽扇豆 RAD测序有效数据统计
2 )将两个个体中相同的短序列进行统计计数， 并过滤掉测序深度为 1的测序数据， 结果统计如表 2所示：
Figure imgf000013_0002
表 2、 羽扇豆 RAD-tag统计 图 10示出 RAD-tag测序深度分布图。
3 )将两个个体计数后的短序列数据内部进行比对， 聚类， 过滤重复 区域的序列， 并生成杂合 RAD-tag的一致性序列。
4 )接着两个个体之间的数据使用比对软件进行比对寻找 INDEL 位 点。 过滤时， 比对容许的错配数为 5, 空位数为 2， 即一个 RAD-tag上最 多容许存在 2个 INDEL位点。
综上， 通过以上步骤的处理， 在羽扇豆父本和母本两个个体中， 总 共找到了 753个 INDEL位点标记。
图 11 示出本发明的基因组 INDEL位点标记装置的一个实施例的结 构图。 如图 11 所示， 该装置包括： 序列深度获取单元 113, 用于统计两 个个体基因组的 RAD 单端测序序列的测序深度； 序列深度过滤单元 115， 用于过滤掉测序深度为 1 的测序序列； 一致性序列形成单元 117, 确定两个个体基因组内部 RAD-tag 上的杂合位点并形成一致性序列。 例 如， 一致性序列形成单元 117对个体基因组内部的 RAG-tag进行不容许 空隙的两两比对以确定个体基因组内部 RAD-tag 上的杂合位点并形成一 致性序列。 不容许空隙的两两比对的容许的错配数根据测序序列的长度 确定。 位点确定单元 119, 用于将两个个体的非杂合 RAD-tag 的序列和 杂合 RAD-tag的一致性序列进行比对， 以确定两个个体在 RAD-tag上的 INDEL信息。
图 12示出本发明的基因组 INDEL位点标记装置的另一个实施例的 结构图。 如图 12所示， 在一个实施例中， 该 INDEL位点标记装置包括 序列深度获取单元 113、 序列深度过滤单元 115、 一致性序列形成单元 117、 重复序列过滤单元 128、 位点确定单元 129。 其中， 序列深度获取单 元 113、 序列深度过滤单元 115、 和一致性序列形成单元 117可以参见上 述实施例的描述， 为简洁起见在此不再详细描述。 重复序列过滤单元 128 用于过滤掉两个个体基因组内部处于重复区域的 RAD-tag。 例如， 满足 如下条件作为处于基因组序列的重复区域中的 RAD-tag序列： 测序序列 在基因组上存在两个以上拷贝， 同时所述两个以上拷贝与对应的同源染 色体在不同的位置上存在杂合位点； 和 /或测序序列在基因组在存在多个 拷贝， 且具有较高的测序深度， 其中一个拷贝与对应的同源染色体上存 在杂合位点。 位点确定单元 129 将两个个体处于非重复区域的非杂合 RAD-tag序列和杂合 RAD-tag的一致性序列进行比对。
在一个实施例中， INDEL 位点标记装置还包括： 测序单元 121， 用 于获得两个个体基因组的 RAD单端测序序列； 测序序列过滤单元 122, 用于对获得的两个个体基因组的 RAD单端测序序列进行过滤以去除不合 格的测序序列。 不合格的测序序列例如包括： 测序质量低于预定的低质 量阈值的碱基个数超过整条测序序列碱基个数的 50%的测序序列； 和 /或 测序序列中测序结果不确定的碱基个数超过整条测序序列碱基个数的 10%的测序序列； 和 /或存在外源序列的测序序列； 和 /或起始的几个碱基 不是酶切末端序列的测序序列。
在一个实施例中， 该装置还包括 INDEL信息过滤单元 130, 对两个个 体基因组在 RAD-tag上的 INDEL信息进行过滤。 例如， INDEL信息过 滤单元 130提取两个个体基因组之间具有 INDEL信息的比对结果； 过滤 掉两个个体基因组之间 RAD-tag上 INDEL数大于 2以及空位长度大于 3 的比对结果； 对每一对的 INDEL 的比对结果， 如果存在其他的比对结 果， 满足比对的错配数小于等于 3并且空位数小于等于 1 , 则将该比对结 果过滤掉； 当 INDEL信息与酶切位点以 SJ^列末端的距离在 3bp以内。 则将该类型的比对结果过滤掉。
图 11、 12 中的单元的功能和实现还可以参见本文中关于方法的描 述。
[参考文献]
1. Li, R. et al. SNP detection for massively parallel whole-genome resequencing. Genome Research 19, 1124 (2009).
2. Li Η·*, Handsaker Β·*, Wysoker A., Fennell T., Ruan J ., Homer N., Marth G" Abecasis G" Durbin R. and 1000 Genome Project Data Processing Subgroup (2009) The Sequence alignment/map (SAM) format and SAMtools. Bioinformatics, 25, 2078-9. [PMID: 19505943]
3. Hohenlohe, P.A. et al. Population genomics of parallel adaptation in threespine stickleback using sequenced RAD tags. PLoS Genet 6, el000862. 至此， 已经详细描述了根据本发明的方法和装置。 为了避免遮蔽本 发明的构思， 没有描述本领域所公知的一些细节。 本领域技术人员根据 上面的描述， 完全可以明白如何实施这里公开的技术方案。
虽然已经通过示例对本发明的一些特定实施例进行了详细说明， 但 是本领域的技术人员应该理解， 以上示例仅是为了进行说明， 而不是为 了限制本发明的范围。 本领域的技术人员应该理解， 可在不脱离本发明 的范围和精神的情况下， 对以上实施例进行修改。 本发明的范围由所附 权利要求来限定。

Claims

权 利 要 求
1. 一种基因组单***缺失（INDEL )位点的标记方法， 其特征在于， 包括：
统计两个个体基因组的限制性内切位点相关 DNA (RAD)单端测序序列 的测序深度；
过滤掉两个个体基因组中测序深度为 1的 RAD单端测序序列； 确定两个个体基因组内部限制性内切位点相关 DNA标记 (RAD-tag)上 的杂合位点并形成一致性序列；
将两个个体基因组的非杂合 RAD-tag序列及杂合 RAD-tag序列的一 致性序列进行比对， 以确定两个个体在 RAD-tag上的 INDEL信息。
2. 根据权利要求 1所述的方法， 其特征在于，
在将两个个体的非杂合 RAD-tag 的序列及杂合 RAD-tag 的一致性序 列进行比对之前包括：
过滤掉两个个体基因组内部处于重复区域的 RAD-tag序列；
所述将两个个体的非杂合 RAD-tag序列及杂合 RAD-tag序列的一致 性序列进行比对包括：
将两个个体处于非重复区域的非杂合 RAD-tag序列及杂合 RAD-tag 序列的一致性序列进行比对。
3. 根据权利要求 2所述的方法， 其特征在于， 满足如下条件作为处于 基因组序列的重复区域中的 RAD-tag序列：
测序序列在基因组上存在两个以上拷贝， 同时所述两个以上拷贝与对 应的同源染色体在不同的位置上存在杂合位点； 和 /或
测序序列在基因组在存在多个拷贝， 且具有较高的测序深度， 其中一 个拷贝与对应的同源染色体上存在杂合位点。
4. 根据权利要求 1至 3中任意一项所述的方法， 其特征在于， 在统计 两个个体基因组的 RAD单端测序序列的测序深度之前包括：
通过高通量测序获得两个个体基因组的 RAD单端测序序列；
对获得的两个个体基因组的 RAD单端测序序列进行过滤以去除不合 格的测序序列。
5. 根据权利要求 4所述的方法， 其特征在于， 所述不合格的测序序列 包括：
测序质量低于预定的低质量阈值的碱基个数超过整条测序序列碱基个 数的 50%的测序序列； 和 /或
测序序列中测序结果不确定的碱基个数超过整条测序序列碱基个数的
10%的测序序列； 和 /或
存在外源序列的测序序列； 和 /或
起始的几个碱基不是酶切末端序列的测序序列。
6. 根据权利要求 1所述的方法， 其特征在于， 所述确定两个个体基因 组内部 RAD-tag上的杂合位点包括：
个体基因组内部的 RAG-tag进行不容许空隙的两两比对以确定个体基 因组内部 RAD-tag上的杂合位点。
7. 根据权利要求 6所述的方法， 其特征在于， 根据测序序列的长度确 定所述不容许空隙的两两比对的容许的错配数。
8. 根据权利要求 6 所述的方法， 其特征在于， 个体基因组内部的 RAG-tag进行不容许空隙的两两比对以确定个体基因组内部 RAD-tag上 的杂合位点包括：
在个体基因组内部的测序序列之间进行不容许空隙的两两比对； 将所有满足比对条件的测序序列进行聚类；
挑选出聚类结果中只有两种测序序列的聚类结果， 该测序序列的位置 即存在杂合位点。
9. 根据权利要求 1所述的方法， 还包括： 对两个个体基因组在 RAD- tag上的 INDEL信息进行过滤。
10. 根据权利要求 9 所述的方法， 其特征在于， 所¾|~两个个体基因 组在 RAD-tag上的 INDEL信息进行过滤包括：
提取两个个体基因组之间具有 INDEL信息的比对结果；
过滤掉两个个体基因组之间 RAD-tag上 INDEL数大于 2以及空位长 度大于 3的比对结果；
对每一对的 INDEL 的比对结果， 如果存在其他的比对结果， 满足比 对的错配数小于等于 3并且空位数小于等于 1 , 则将该比对结果过滤掉； 当 INDEL信息与酶切位点以 ^列末端的距离在 3bp 以内。 则将该 类型的比对结果过滤掉。
11. 一种基因组单***缺失 （INDEL )位点的标记装置， 其特征在 于， 包括：
序列深度获取单元， 用于统计两个个体基因组的限制性内切位点相关 DNA(RAD)单端测序序列的测序深度；
序列深度过滤单元， 用于过滤掉两个个体基因组中测序深度为 1 的 RAD单端测序序列；
一致性序列形成单元， 用于确定两个个体基因组内部限制性内切位点 相关 DNA标记 (RAD-tag)上的杂合位点并形成一致性序列；
位点确定单元， 用于将两个个体的非杂合 RAD-tag 的序列及杂合 RAD-tag的一致性序列进行比对以确定两个个体在 RAD-tag上的 INDEL 息
12. 根据权利要求 11所述的装置， 其特征在于， 还包括：
重复序列过滤单元， 用于过滤掉两个个体基因组内部处于重复区域的
RAD-tag;
所述位点确定单元将两个个体处于非重复区域的非杂合 RAD-tag序列 及杂合 RAD-tag的一致性序列进行比对以确定两个个体在 RAD-tag上的 INDEL信息。
13. 根据权利要求 12所述的装置， 其特征在于， 满足如下条件作为处 于基因组序列的重复区域中的 RAD-tag序列：
测序序列在基因组上存在两个以上拷贝， 同时所述两个以上拷贝与对 应的同源染色体在不同的位置上存在杂合位点； 和 /或
测序序列在基因组在存在多个拷贝， 且具有较高的测序深度， 其中一 个拷贝与对应的同源染色体上存在杂合位点。
14. 根据权利要求 11至 13 中任意一项所述的装置， 其特征在于， 还 包括：
测序单元， 用于通过高通量测序获得两个个体基因组的 RAD单端测 序序列；
测序序列过滤单元， 用于对获得的两个个体基因组的 RAD单端测序 序列进行过滤以去除不合格的测序序列。
15. 根据权利要求 14所述的装置， 其特征在于， 所述不合格的测序序 列包括：
测序质量低于预定的低质量阈值的碱基个数超过整条测序序列碱基个 数的 50%的测序序列； 和 /或
测序序列中测序结果不确定的碱基个数超过整条测序序列碱基个数的 10%的测序序列； 和 /或
存在外源序列的测序序列； 和 /或
起始的几个碱基不是酶切末端序列的测序序列。
16. 根据权利要求 11所述的装置， 其特征在于， 一致性序列形成单元 对个体基因组内部的 RAG-tag 进行不容许空隙的两两比对以确定个体基 因组内部 RAD-tag上的杂合位点并形成一致性序列。
17. 根据权利要求 16所述的装置， 其特征在于， 所述不容许空隙的两 两比对的容许的错配数根据测序序列的长度确定。
18. 根据权利要求 11所述的装置， 还包括：
INDEL 信息过滤单元， 用于对两个个体基因组在 RAD-tag 上的 INDEL信息进行过滤。
19.根据权利要求 18所述的装置， 其特征在于， 所述 INDEL信息过 滤单元提取两个个体基因组之间具有 INDEL信息的比对结果； 过滤掉两 个个体基因组之间 RAD-tag上 INDEL数大于 2以及空位长度大于 3的比 对结果； 对每一对的 INDEL 的比对结果， 如果存在其他的比对结果， 满 足比对的错配数小于等于 3并且空位数小于等于 1， 则将该比对结果过滤 掉； 当 INDEL信息与酶切位点以 Sjf列末端的距离在 3bp以内。 则将该 类型的比对结果过滤掉。
PCT/CN2012/071329 2011-12-29 2012-02-20 基因组indel位点标记方法和装置 WO2013097328A1 (zh)

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
CN201110448880.5 2011-12-29
CN201110448880 2011-12-29

Publications (1)

Publication Number Publication Date
WO2013097328A1 true WO2013097328A1 (zh) 2013-07-04

Family

ID=48696314

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/CN2012/071329 WO2013097328A1 (zh) 2011-12-29 2012-02-20 基因组indel位点标记方法和装置

Country Status (1)

Country Link
WO (1) WO2013097328A1 (zh)

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105779433A (zh) * 2014-12-15 2016-07-20 天津华大基因科技有限公司 试剂盒及其用途
CN105779434A (zh) * 2014-12-15 2016-07-20 天津华大基因科技有限公司 试剂盒及其用途
CN105779435A (zh) * 2014-12-15 2016-07-20 天津华大基因科技有限公司 试剂盒及其用途
CN105779432A (zh) * 2014-12-15 2016-07-20 天津华大基因科技有限公司 试剂盒及其用途

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
BAIRD, N.A. ET AL.: "Rapid SNP Discovery and Genetic Mapping Using Sequenced RAD markers", PLOS ONE, vol. 3, no. 10, October 2008 (2008-10-01), pages E3376 *
LI, DONG: "Constructing High-resolution Genetic Variation Maps and Related Studies on Silkworm Nuclear Genomes and Silkworm Mitochondrial Genomes", AGRICULTURE, CHINA DOCTORAL DISSERTATIONS FULL-TEXT DATABASE, 15 September 2011 (2011-09-15), pages 19 - 21 *
WF, P. ET AL.: "Mapping with RAD (restriction-site associated DNA) markers to rapidly identify QTL for stem rust resistance in Lolium perenne", TAG THEORETICAL AND APPLIED GENETICS, vol. 122, no. 8, 23 February 2011 (2011-02-23), pages 1467 - 1480 *

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105779433A (zh) * 2014-12-15 2016-07-20 天津华大基因科技有限公司 试剂盒及其用途
CN105779434A (zh) * 2014-12-15 2016-07-20 天津华大基因科技有限公司 试剂盒及其用途
CN105779435A (zh) * 2014-12-15 2016-07-20 天津华大基因科技有限公司 试剂盒及其用途
CN105779432A (zh) * 2014-12-15 2016-07-20 天津华大基因科技有限公司 试剂盒及其用途

Similar Documents

Publication Publication Date Title
CN108350494B (zh) 用于基因组分析的***和方法
CN110997937B (zh) 具有可变长度非随机独特分子标识符的通用短衔接子
EP3271480B1 (en) Screening for structural variants
JP2019523638A (ja) 遺伝子突然変異を検出するマルチポジショニングダブルタグアダプターセット、及びその調製方法と応用
WO2013097048A1 (zh) 基因组单核苷酸多态性位点的标记方法和装置
CN111755072B (zh) 一种同时检测甲基化水平、基因组变异和***片段的方法及装置
WO2014197377A2 (en) Methods and systems for storing sequence read data
EP2834762A1 (en) Sequence assembly
CN106715711A (zh) 确定探针序列的方法和基因组结构变异的检测方法
US20210375397A1 (en) Methods and systems for determining fusion events
Larson et al. A clinician’s guide to bioinformatics for next-generation sequencing
Guo et al. Single-nucleotide variants in human RNA: RNA editing and beyond
WO2013097328A1 (zh) 基因组indel位点标记方法和装置
JP2023523002A (ja) 染色体近接実験における構造的変異検出
WO2012097474A1 (zh) 检测转基因外源片段***位点的方法和***
Jiang et al. Long-read based novel sequence insertion detection with rCANID
US20200395098A1 (en) Alignment using homopolymer-collapsed sequencing reads
CN114420213A (zh) 一种生物信息分析方法及装置、电子设备及存储介质
JP2008161056A (ja) Ｄｎａ配列解析装置、ｄｎａ配列解析方法およびプログラム
CN102831331A (zh) 基于酶切建库双末端测序的长度多态性标记的引物设计开发方法
KR20220064959A (ko) 낮은 빈도 변이의 검출 및 리포팅을 용이하게 하기 위한 dna 라이브러리 생성 방법
US20210164033A1 (en) Method and system for nucleic acid sequencing
KR101977976B1 (ko) 앰플리콘 기반 차세대 염기서열 분석기법에서 프라이머 서열을 제거하여 분석의 정확도를 높이는 방법
Ning et al. ssahaSNP-a polymorphism detection tool on a whole genome scale
WO2013097143A1 (zh) 估计基因组杂合率的方法和装置

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 12861668

Country of ref document: EP

Kind code of ref document: A1

NENP Non-entry into the national phase

Ref country code: DE

32PN Ep: public notification in the ep bulletin as address of the adressee cannot be established

Free format text: NOTING OF LOSS OF RIGHTS PURSUANT TO RULE 112(1) EPC (EPO FORM 1205A DATED 07/11/2014)

122 Ep: pct application non-entry in european phase

Ref document number: 12861668

Country of ref document: EP

Kind code of ref document: A1