A method for splicing first-generation sequencing data for 16S rDNA identification

Abstract

The present invention discloses a method for splicing first-generation sequencing data for 16S rDNA identification. It includes: establishing a sequence splicing program based on an alignment module and a splicing module, realizing the splicing of sequence data to be analyzed under the first-generation sequencing machine, and outputting a complete sequence that can be used for subsequent identification and analysis; the alignment module includes an alignment model, a scoring standard and a preferred standard; the alignment model uses one direction for scoring calculation and backtracking; the scoring standard introduces the concept of degenerate bases, and expands 4 bases to 15 bases, including A, T, C, G, R, Y, M, K, S, W, H, B, V, D, N; wherein, when the completely identical bases and degenerate bases contain normal bases, the score is 10 points, any base alignment N is 1 point, any base alignment "-" is 2 points, and two bases in the same site are different and are 2 points. The present invention establishes the core idea of ​​brand-new sequence splicing, with the advantages of fast speed, small amount of calculation, and high accuracy.

Description

Technical Field

[0001] The present invention relates to the field of microbial technology, and in particular to a sequence splicing algorithm for first-generation sequencing data in a 16S rDNA sequencing experiment step. Background Art

[0002] The splicing of short sequences is a necessary condition for completing the "human microbiome engineering" and is the key to the promotion and application of the next-generation sequencing technology. Because if you want to obtain the whole genome or partial genome sequence of the microorganisms in the sample, you must use sequence splicing technology to complete it. Based on the current sequencing technology, the short sequences obtained are much smaller than the original genome size and cannot provide directly usable knowledge to the relevant scientists. Therefore, the obtained data needs to be spliced. The so-called splicing is to design an algorithm to form overlapping groups (contigs), long sequences with gaps (scaffolds) or even entire genomes with the obtained short sequences. Sequence splicing is the basis of genetic information research and the first processing step in a large number of sequence analysis processes. Its accuracy is crucial for all downstream analyses.

[0003] During the implementation of the Human Genome Project, sequence splicing algorithms were born one after another. The general idea of ​​the algorithm based on the Greedy-extension strategy is as follows: starting from the short sequence with the largest overlap, the short sequences that meet the overlap conditions are spliced ​​into contigs until there are no reads or contigs to continue splicing. The greedy algorithm uses the local optimal overlap to splice reads locally. The optimal judgment criteria vary for different algorithms, but most of them are measured by overlap length and similarity. The algorithm strategy based on Overlap-Layout-Consensus is different from the greedy algorithm. It considers the overlap part globally and analyzes all reads to find the optimal splicing. The first step is the same as the greedy strategy, calculating the overlap between reads. The overlap graph is constructed based on the overlap information. Each read represents a node. If two reads overlap, the two reads are connected by an edge. Finally, the graph is analyzed to find the splicing path of the reads. The important step of the Overlap-Layout-Consensus algorithm is the last step, namely the Layout stage, which is the core of the OLC strategy and determines the final splicing result. The algorithm based on de Bruijn graph is also called Euler path or kmer graph algorithm. The idea of ​​de Bruijn graph splicing is inspired by hybrid sequencing experiments. The sequencing results of the experiment are a series of fixed-length gene sequence fragments from the genome. These fragments are kmers. These fragments are processed and then spliced. Kmer is obtained by sliding segmentation of reads, starting from each position of the read, and sliding segmentation is divided into fragments of k bases. Summary of the invention

[0004] In order to overcome the shortcomings of the prior art, the purpose of the present invention is to provide a first-generation sequencing data splicing method for 16S rDNA identification, which is applied to the splicing process of 16S rDNA sequencing analysis, can significantly improve the sequence splicing efficiency and reduce the amount of program calculations.

[0005] Specifically, the present invention provides the following technical solutions:

[0006] The present invention provides a sequence splicing method applied to a 16S rDNA sequencing analysis process. The input data is the downstream sequence data of the first-generation sequencing. Based on the software of the present invention, a complete sequence for downstream analysis can be obtained quickly and accurately.

[0007] The gene sequence splicing software provided by the present invention is based on the core idea of ​​global comparison, and further changes and optimizes on this basis, uniquely adds degenerate bases to the operation, and realizes the purpose of sequence comparison and splicing according to specific scoring standards and operation rules. The specific core steps are as follows: the first step is to input the primer and the original machine sequence into the comparison module to complete the comparison, and according to the scoring standard of the present invention, sort according to specific indicators, and select the best score sequence as the seed sequence; the second step is to input the seed sequence and other original machine sequences into the comparison module to complete the comparison, and according to the same scoring standard, sort according to specific indicators to select the best, and then input the seed sequence and the best sequence of this comparison into the splicing module to realize the first splicing; the third step is to use the result of the first splicing as the seed sequence, and repeat the second step to realize the comparison and selection; the fourth step is to splice the seed sequence with the preferred sequence, and if there are more sequences that need to be spliced, the number of splicing is counted based on the primer information, and then the recursive splicing mode is turned into according to the second step until the splicing is completed.

[0008] Regarding the scoring criteria, the present invention specifically introduces the concept of degenerate bases, expands 4 bases to 15 bases, and performs matching scoring within a certain range of rules (when degenerate bases contain different bases). Specifically, if two sequences being compared have exactly the same normal bases at the same position, the score is 10 points; if the degenerate bases contain normal bases, the score is 10 points (such as R bases are compared to normal bases A / G, the score is 10 points); N is compared to any base for 1 point; any base is compared to "-" for 2 points; two bases at the same position are completely different for 2 points.

[0009] Specifically, the comparison module of the present invention implements the following steps: (1) Assume that the two input sequences (which can be regarded as two strings with fixed characters) are s and t, and its main data structure is the s*t matrix a. After initialization according to the scoring standard of the present invention, the matrix a will obtain new values ​​in turn. Here, two of the three directions of the traditional local comparison are abandoned, and only one direction is followed to calculate according to the following formula (1):

[0010] (1);

[0011] Where a represents the matrix, a[i,j] is the score of the current site, g is the gap penalty, a[i-1,ji] represents the score of the upper left corner, p(i,j) represents the current score, j and i are s[1...i] and t[1....j] respectively. (2) This module selects the maximum value in the last row or last column of matrix a as the backtracking starting point and backtracks in a unique direction. The backtracking calculation formula (2) is as follows:

[0012] (2);

[0013] Where Sim(s,t) represents the maximum score. After calculating the maximum score, use it as the starting point for backtracking and follow the principles of the dynamic programming algorithm to backtrack to a[0,0].

[0014] The global alignment method used in the present invention is different from the traditional alignment method based on the above-mentioned scoring criteria. Traditional global alignment or semi-global alignment cannot avoid gaps in the middle of the sequence. In order to avoid the gaps commonly used in sequence alignment, the traditional scoring method is changed and combined with the algorithm idea of ​​dynamic programming to fit the optimal seed sequence.

[0015] Among them, regarding the selection criteria for sequences during sequence alignment, the present invention calculates scores and sorts according to the two indicators of "including degenerate bases and base consistency" and "strict base consistency" in the alignment region, and the two indicators are X and Y. When sorting, the present invention gives priority to the score of the "including degenerate bases and base consistency" indicator. If the scores of multiple sequences are the same, they are randomly selected. The specific calculation formula is as follows:

[0016] (3);

[0017] Specifically, the splicing module of the present invention is implemented as follows: the two bases corresponding to each site of the two sequences are compared, and the comparison result is calculated using the following formula:

[0018] (4);

[0019] Among them, f(a,b) and f(a,b,d) represent degenerate bases. f(x) represents the output bases, f(a) and f(b) are common ATCG bases, and '-' represents the gap filler.

[0020] In addition, the present invention follows the following rules during the splicing process: 1) When any base corresponds to the N base or the '-' symbol, a more specific base is selected; 2) For two degenerate bases containing the same base, the intersecting base or the degenerate base is selected; 3) When two degenerate bases do not have intersecting bases, the degenerate bases of the union are taken; 4) When the 'A', 'T', 'C', and 'G' bases are the same, the corresponding bases are taken; if they are different, the degenerate bases of the union are taken.

[0021] The beneficial effects of the present invention are:

[0022] The sequence splicing software of the present invention has excellent performance in sensitivity, accuracy, specificity and precision. For similar splicing software, too high or too low sensitivity will cause accuracy problems, and the relationship between the two is like a normal distribution. If the sensitivity is too high, gaps are easily generated during the comparison process, and most bases cannot be successfully compared to produce the optimal solution, so the accuracy will be affected. And if the sensitivity is too low, the random probability of the comparison will also be increased, thereby reducing the accuracy of the sequence comparison. Therefore, the sequence splicing algorithm of the present invention takes this into account when using the scoring matrix and the penalty matrix, and the penalty mechanism is improved. In a slightly modified manner, the penalty is changed to a cumulative score for a specific situation, and the maximum value of the scoring matrix is ​​obtained, so that the problem of increasing gaps and increasing random probability caused by too many penalties can be avoided, thereby reducing errors. In addition, the present invention adopts an improved dynamic programming algorithm, using degenerate bases as the final splicing method, and compared with other software that does not support degenerate bases, the result produced by the present invention has more accurate characteristics when describing bases. BRIEF DESCRIPTION OF THE DRAWINGS

[0023] Figure 1 Shown is a schematic flow chart of a method for splicing 16S rDNA sequencing first generation data in an embodiment of the present invention. DETAILED DESCRIPTION

[0024] The following examples are used to illustrate the present invention but are not intended to limit the scope of the present invention.

[0025] Example 1: Splicing of 16S rDNA Sequencing Data

[0026] This embodiment provides a sequence splicing algorithm applied to the 16S rDNA sequencing analysis process, and uses the Python language to implement the program based on the program. Based on the program, 2253 test samples from the first-generation sequencing of 16S sequencing analysis were spliced ​​to obtain a complete sequence for subsequent species identification. Figure 1 As shown, the specific implementation steps are as follows:

[0027] (1) The data from the first-generation sequencing technology include ab1 files and seq files. The input data required by the implementation program of the present invention is a fasta file composed of seq data. The specific parameters are as follows: --barcode Primer used for sequencing, use the default "GTGCCAGCAGCCGCGGTAA"; --output Output file path; --process F / T Whether to output the detailed program operation process, set to F.

[0028] (2) Based on the fasta file and barcode information prepared in step (1), run the program to perform sequence assembly. Obtain a complete assembled fasta sequence file assembly.fasta, which can be used for subsequent species identification analysis.

[0029] (3) Use phrap software to splice 2253 test samples. The specific key parameters are: -minmatch 20.

[0030] (4) Use clustalw and needle software to identify species of the spliced ​​sequences obtained by the above two splicing methods. The default parameters of clustalw were set. The key parameters of needle software are: -aformat3 fasta -gapopen10.0 -gapextend 0.5. The final species identification results were obtained based on the similarity sorting.

[0031] The sequence assembly and species identification results of the test samples are shown in the following table:

[0032]

[0033] The splicing results show that this method can accurately and quickly realize the splicing of large quantities of 16S rDNA first-generation sequencing data. Compared with other methods, this method has 87.8% of the samples that are exactly the same as Phrap software in terms of splicing length. For the remaining 12.2% of samples with different splicing results, there are still 194 samples that can be consistent with Phrap in species identification results.

[0034] Although the present invention has been described in detail above with general descriptions and specific embodiments, it is obvious to those skilled in the art that some modifications or improvements may be made thereto based on the present invention. Therefore, these modifications or improvements made without departing from the spirit of the present invention all fall within the scope of protection claimed by the present invention.

Claims

1. A method for splicing first generation sequencing data for 16S rDNA identification, characterized in that: include: Establish a sequence splicing program based on the alignment module and the splicing module to realize the splicing of the sequence data to be analyzed after the first-generation sequencing machine, and output a complete sequence that can be used for subsequent identification and analysis; The comparison module includes a comparison model, a scoring standard and a selection standard; The comparison model uses one direction for scoring calculation and backtracking; The scoring standard introduces the concept of degenerate bases, expanding the 4 bases to 15 bases, including A, T, C, G, R, Y, M, K, S, W, H, B, V, D, and N. Among them, the score is 10 points when the bases are completely identical and the degenerate bases contain normal bases, 1 point when any base pair is N, 2 points when any base pair is "-", and 2 points when two bases at the same site are different; The optimization criteria are to calculate the scores and sort according to the two indicators of "including degenerate bases and base consistency" and "strict consistency of bases" in the comparison region, and the score of "including degenerate bases and base consistency" is given priority in sorting; the calculation formula (3) is as follows: （3）； In the i-th base alignment, X is the statistical score of the "including degenerate bases and base consistency" indicator, and Y is the statistical score of the "strict consistency of bases" indicator; The splicing module includes a splicing algorithm model and splicing rules followed by the model; The splicing algorithm model adds degenerate bases for operation, and the output splicing result uses degenerate bases to represent uncertain bases in the overlapping region; The splicing rules followed by the splicing algorithm model include: when any base corresponds to the N or '-' symbol, select a more specific base; when two degenerate bases contain the same base, select the intersecting base or the degenerate base; when two degenerate bases do not have intersecting bases, take the degenerate base of the union; when the 'A', 'T', 'C', 'G' bases are the same, take the corresponding base, and if they are different, take the degenerate base of the union; In the first step, the primers and the original sequence are input into the comparison module to complete the comparison, and the comparison area is sorted according to the two indicators of "including degenerate bases and base consistency" and "strict consistency of bases" according to the scoring standard, and the sequence with the best score is selected as the seed sequence; in the second step, the seed sequence and other original sequences are input into the comparison module to complete the comparison, and the selection is made according to the specific indicators according to the same scoring standard, and then the seed sequence and the best sequence of this comparison are input into the splicing module to realize the first splicing; in the third step, the result of the first splicing is used as the seed sequence, and the second step is repeated to realize the comparison and selection; in the fourth step, the seed sequence and the preferred sequence are spliced. If there are more sequences that need to be spliced, the number of splicing is counted based on the primer information, and then the recursive splicing mode is switched to according to the second step until the splicing is completed; The comparison module performs the following steps: 1) Assume that the two input sequences are s and t, and the data structure is a matrix a of s*t; after being initialized according to the scoring criteria, the matrix a will obtain new values ​​in turn and be calculated according to the following formula (1): （1）； Where a represents the matrix, a[i,j] is the score of the current site, g is the gap penalty, a[i-1,ji] represents the score of the upper left corner, p(i,j) represents the current score, j and i are s[1...i] and t[1....j] respectively; 2) This module selects the maximum value in the last row or last column of matrix a as the backtracking starting point and backtracks in a unique direction. The backtracking calculation formula (2) is as follows: （2）； Among them, Sim(s,t) represents the maximum score. After calculating the maximum score, it is used as the backtracking starting point and backtracking is performed according to the principles of the dynamic programming algorithm until a[0,0]; The splicing module is implemented as follows: the two bases corresponding to each site of the two sequences are compared, and the comparison result is calculated using the following formula (4): （4）； Among them, f(a,b) and f(a,b,d) represent degenerate bases; f(x) represents the output base, f(a) and f(b) are common ATCG bases, and '-' represents the gap filler. During the splicing process, the following rules are followed: a) When any base corresponds to the N base or the '-' symbol, the more specific base is selected; b) For two degenerate bases that contain the same base, the intersecting base or the degenerate base is selected; c) When two degenerate bases do not have intersecting bases, the degenerate bases of the union are taken; d) When the 'A', 'T', 'C', and 'G' bases are the same, the corresponding bases are taken; if they are different, the degenerate bases of the union are taken.

Images