Nucleic acid sequence clustering method, apparatus, computer readable storage medium, terminal

Abstract

The application discloses a nucleic acid sequence clustering method and device, a computer readable storage medium and a terminal. The terminal constructs a tree structure with multiple branches to search a specified interval of a nucleic acid sequence, thereby avoiding a large amount of time consumed by traditional calculation of an editing distance. In addition, the application adopts a node drift algorithm to resist interference caused by errors in the nucleic acid sequence. Compared with existing nucleic acid clustering algorithms, the method provided by the application can cluster a large number of unidentified nucleic acid sequences, and has the functions of automatically correcting and comparing the clustered nucleic acid sequences. The method can directly output the corrected nucleic acid original sequence, thereby greatly reducing the processing time after sequencing reading.

Description

TECHNICAL FIELD

[0001] The present application belongs to the technical field of data storage, and particularly relates to a nucleic acid sequence clustering method, device, computer readable storage medium and terminal. BACKGROUND

[0002] Nucleic acid is the general term of deoxyribonucleic acid (DNA) and ribonucleic acid (RNA), which is a biological macromolecular compound polymerized by many nucleotide monomers and is one of the most basic substances of life. The research of nucleic acid involves multiple fields such as biology, medicine and computer.

[0003] Sequencing is one of the basic means for studying nucleic acid. Through sequencing technology, the information of nucleic acid molecules can be read into a computer storage medium for further use and analysis. In recent years, with the maturity of the second-generation sequencing technology, the related field research has entered a high-speed development. The second-generation sequencing (NGS) is also called high-throughput sequencing, which initiatively introduces a reversible terminal to realize sequencing by synthesis, and determines the nucleic acid sequence by capturing the special markers (primers) carried by the newly added bases in the nucleic acid replication process. The second-generation sequencing has two important features: 1. high throughput, the second-generation sequencing can parallelly sequence several million nucleic acid molecules at a time; 2. short sequence length, since the sequencing process decreases with the growth of read length, the cooperativity of gene cluster replication decreases, which leads to the decline of sequencing quality, so the read length of the second-generation sequencing is not more than 500 bp. For longer genomes and metagenomes, they need to be broken into small fragments for sequencing, and then spliced after sequencing.

[0004] Due to the large number of nucleic acid sequences generated during the sequencing of the second-generation sequencing technology, it will be difficult to analyze the next step, so it is necessary to use a clustering algorithm to classify the sequenced sequences to obtain the original sequence cluster for analysis to improve the efficiency, Figure 2 The flowchart of nucleic acid sequence sequencing reading is shown. However, in the field of DNA storage and the like, in order to improve the accuracy of reading sequences, multiple rounds of molecular chain amplification and deep sequencing are often performed, which can produce DNA sequences up to hundreds of millions. The existing clustering method will take an intolerable time and memory occupation. Therefore, the clustering analysis method for a large number of nucleic acid sequences still needs to be improved. CN 110111843 A provides a method for clustering nucleic acid, but since it uses the method of calculating the edit distance between sequences to cluster, a large amount of calculation time is still required, which is difficult to process quickly for complex sequencing data.

[0005] The difficulties of the current nucleic acid clustering algorithm mainly concentrate on the following aspects:

[0006] (1) The number of nucleic acid sequences that need to be clustered is huge, often with tens of millions or even higher orders of magnitude, and the number of clusters is extremely large. The traditional clustering algorithm needs to consume a lot of time and memory. At present, when the traditional clustering algorithm is used for more than 100,000 sequences, the time required will exceed 10 hours, and the time required will increase sharply with the increase of the number of sequences. It is almost impossible to process DNA sequences of several million.

[0007] (2) As a text form of sequence string, nucleic acid sequence cannot use the traditional clustering algorithm based on Euclidean distance. The existing nucleic acid clustering algorithm still uses edit distance to determine the distance, so it has high computational complexity. There is still no nucleic acid clustering algorithm with linear computational complexity.

[0008] (3) Nucleic acid sequences will randomly produce errors during amplification and sequencing, which will inevitably increase the difficulty of clustering. For some high error rate sequences, the existing algorithm is difficult to classify correctly. SUMMARY

[0009] In view of the technical problems existing in the prior art, the purpose of the present application is to provide a nucleic acid sequence clustering method, device, computer readable storage medium and terminal, which constructs a tree structure with multiple branches to search the specified interval of the nucleic acid sequence, thereby avoiding the large amount of time consumed by traditional calculation of edit distance.

[0010] To achieve the purpose of the present application, the technical scheme provided by the present application is as follows:

[0011] First aspect

[0012] The present application provides a nucleic acid sequence clustering method, comprising the following steps:

[0013] Step a: initialize parameters;

[0014] Step c: determine whether the test sequence is empty, if yes, go to step d; if no, go to step f;

[0015] Step d: output the cluster situation and the core sequence set;

[0016] Step e: end;

[0017] Step f: retrieve a test sequence from the tree structure;

[0018] Step g: determine whether a similar core sequence can be retrieved, if yes, go to step h; if no, go to step j;

[0019] Step h: classify the test sequence and the similar core sequence into the same cluster;

[0020] Step i: optionally perform global alignment and correct core sequence errors;

[0021] Step j: add the sequence to be sequenced to the core sequence set and the tree structure, and jump to step c.

[0022] In step a, the initialization parameters include the number of tree structures, the length of the tree structure, the interval position selected by the tree structure, the vertical drift value, the horizontal drift value threshold, the number of processes, the input file format, the output file format, the core sequence set, and the tree structure.

[0023] If it is a multi-process mode, after step a and before step c, step b is further included: shunting the sequence to be sequenced.

[0024] In steps f-j, the following steps are included:

[0025] A sequence to be sequenced is taken out and searched with the tree structure to determine whether a similar core sequence can be searched. First, whether the first interval and the first tree structure can be matched is compared. If they can be matched, the subsequent tree structure search is stopped, the sequence to be sequenced is classified into the same cluster as the similar core sequence, and if the global alignment function is turned on, the sequence is aligned with the matched sequence at this time. After global alignment, the different base positions in the sequence will be marked, and if a base position of a sequence in the core sequence set is frequently marked, it will be considered as an error base position and corrected. If the first interval cannot be matched, the end interval and the end tree structure are matched. If the end interval can be successfully matched, the subsequent tree structure search is stopped, the sequence to be sequenced is classified into the same cluster as the similar core sequence, and if the global alignment function is turned on, the sequence is aligned with the matched sequence at this time. After global alignment, the different base positions in the sequence will be marked, and if a base position of a sequence in the core sequence set is frequently marked, it will be considered as an error base position and corrected. If the end interval cannot be successfully matched, the middle interval is matched. When the middle interval is selected, the number of base positions allowed to drift horizontally and vertically is allowed to drift based on the original specified interval, and then multiple intervals are allowed to be selected. After searching the base interval tree structure, the interval with the smallest horizontal drift value is selected and matched with the matched sequence. If the horizontal drift value is less than the threshold value at this time, the matching is successful. If the global alignment function is turned on, the sequence is aligned with the matched sequence at this time. After global alignment, the different base positions in the sequence will be marked, and if a base position of a sequence in the core sequence set is frequently marked, it will be considered as an error base position and corrected. If it still cannot be matched, the sequencing sequence is considered as a new cluster, added to the core sequence set, and the specified interval is added to the tree structure.

[0026] Wherein, between step c and step d, further comprising step k, setting a threshold, if a cluster contains less sequences, considering the cluster as a noise cluster, discarding it from the core sequence set, in addition, if selecting multi-process mode, merging the core sequence sets between different processes, but different sequences in the sequence set will not be merged, if the input file is a labeled dataset, calculating the time-consuming, accuracy and correctness, finally outputting the cluster result and the core sequence set.

[0027] Second aspect

[0028] The application provides a nucleic acid sequence clustering device, comprising the following units:

[0029] Parameter initialization unit, used for initializing parameters;

[0030] Shunting unit, used for shunting the to-be-tested sequences;

[0031] To-be-tested sequence judgment unit, used for judging whether the to-be-tested sequence is empty;

[0032] Result output unit, used for outputting the cluster situation and the core sequence set;

[0033] Retrieval unit, used for retrieving a to-be-tested sequence with a tree structure;

[0034] Retrieval result judgment unit, used for judging whether a similar core sequence can be retrieved;

[0035] Judgment result division unit, used for dividing the to-be-tested sequence and the similar core sequence into the same cluster, or adding the to-be-tested sequence into the core sequence set and the tree structure;

[0036] Global alignment error correction unit, used for performing global alignment and correcting core sequence errors.

[0037] Third aspect

[0038] The application provides a computer readable storage medium, comprising programs or instructions, when the programs or instructions are run on a computer, realizing any nucleic acid sequence clustering method as described above.

[0039] Fourth aspect

[0040] The application provides a computer terminal, comprising a memory, and one or more processors in communication connection with the memory;

[0041] The memory stores instructions executable by the one or more processors, and the instructions are executed by the one or more processors, so that the one or more processors realize the nucleic acid sequence clustering method as described above.

[0042] Compared with the prior art, the method of the present application bypasses the conventional way of calculating the edit distance, and classifies the plurality of nucleic acid sequences based on the difference in the arrangement of bases in the specified interval between the plurality of nucleic acid sequences to determine the correct original cluster set. When different nucleic acid sequences are classified by comparison, the method allows the movement between the nodes of the tree structure in the storage medium to resist slight differences in sequences within the same cluster, thereby improving the accuracy of classification. The present application allows the classification of sequences to be performed on all bases to improve the accuracy of sequence classification, while allowing the output of original sequences after classification. On this basis, further provided are a device for clustering nucleic acid sequences and a computer-readable storage medium, and a terminal. The method and device of the present application can quickly classify a large number of nucleic acid sequences and obtain original clusters and original sequences, and then perform subsequent analysis of nucleic acid sequences. BRIEF DESCRIPTION OF DRAWINGS

[0043] Figure 1 Flowchart of the nucleic acid sequence clustering method provided in the present application;

[0044] Figure 2 Flowchart of nucleic acid sequence sequencing reads in the prior art;

[0045] Figure 3 Schematic diagram of the tree structure composed of DNA in the present application. DETAILED DESCRIPTION

[0046] It should be noted that the embodiments in the present application and the features in the embodiments can be combined with each other without conflict.

[0047] The present application will be further described in detail below in conjunction with the drawings and specific embodiments. It should be understood that the specific embodiments described herein are only used to explain the present application and do not limit the present application.

[0048] It needs to be explained that the present application relates to gene sequencing, and the specific problem to be solved is to cluster the existing sequencing data, and then restore the original sequence information, and reduce the difficulty of processing sequencing information. When DNA data is in a molecular state, it needs to be read by a sequencer. However, before reading, the DNA molecular chain needs to be copied many times, and after reading, a large amount of repeated DNA information will be generated. Therefore, a clustering software is needed to classify the same type of DNA information into a class to facilitate information reading. For example, there are 10 DNA chains (molecular state), which need to be amplified and copied 10 times, and then thrown into the sequencer to read the information of the DNA chain. There are 100 DNA chain information in the computer, but many of them are repeated, which affects the use. Therefore, a clustering method is needed to re-cluster them into 10 categories for convenient use. (In practice, there will be tens of millions or even hundreds of millions of chains). The nucleic acid clustering method provided by the present application can cluster a large number of unidentified nucleic acid sequences, and has the functions of automatic error correction and comparison of the clustered nucleic acid sequences, and can directly output the corrected nucleic acid original sequence, thereby greatly reducing the processing time after sequencing reading.

[0049] First, the definition of the sequence to be sequenced in the present application is given, and the sequence to be sequenced represents the nucleic acid sequence that has not been classified after sequencing processing. The model of the present application can be briefly described as follows: first, the model will construct a core sequence set (the core sequence set is empty at the beginning), and then each sequence in the sequence to be sequenced is compared with the core sequence set. If the sequence can be compared with each sequence in the core sequence set, it is successfully classified into a specified cluster, otherwise it is added to the core sequence set as a new core sequence. When comparing the sequence in the sequence to be sequenced with the core sequence set, the core sequence set is first constructed as an index of a tree structure, and then compared with the sequence to be sequenced, thereby avoiding the problem of increasing comparison time caused by the increase of the core sequence set. In addition, node drift is allowed on the tree during retrieval to reduce the influence of sequence errors. The tree structure and node drift algorithm are the key of the present application. The following will first introduce these two points, and then give the whole process and description of the method.

[0050] Tree structure

[0051] The tree structure is an important nonlinear data structure. It is a structure in which data elements are organized in a branching relationship, much like a tree in nature. The definition thereof is given herein:

[0052] A tree is a finite set of n (n>0) elements, wherein:

[0053] (1) each element is called a node;

[0054] (2) There is a specific node, called root node or root;

[0055] (3) In addition to the root node, the remaining nodes are divided into m (m≥0) non-intersecting finite sets, and each subset is a tree.

[0056] Figure 3 The schematic diagram of the tree structure for DNA sequence is shown. Since the DNA sequence can only be composed of {A, T, G, C}, each node under the root has at most 4 nodes. The depth of a tree is defined as L if the tree is constructed by n (n>0) sequences with length L. Therefore, a theorem can be obtained as follows:

[0057] For a tree with depth L constructed by M (m≥0) sequences, the computational complexity of any sequence searching the tree is O(L).

[0058] It can be seen that the time for searching the tree nodes is not affected by the number of sequences (nodes) contained in the tree structure.

[0059] Node drift

[0060] During the amplification and sequencing of nucleic acid molecules, errors such as base loss, addition, substitution, etc. may occur randomly. In order to resist the interference of error bases in the nucleic acid sequence after sequencing, the present application allows a certain degree of node drift when searching the sequence to be sequenced in the tree structure to prevent the error sequence from being successfully matched to the correct core sequence. The drift is divided into horizontal drift and vertical drift:

[0061] Horizontal drift: When a specified node does not exist under a root of the search tree, the remaining existing nodes under the root will be searched. If there are other nodes, and the next node can still be matched when the node is used as a subtree, then drift to the other node is allowed.

[0062] Vertical drift: For a sequence interval [a, b] (b>a≥0) that should be searched in the tree, when the sequence vertically drifts by t (t≤a), then all sliding window intervals between [a+t, b+t] and [a-t, b-t] are used to search the tree.

[0063] As defined, horizontal drift can reduce the search impact caused by base substitution errors. Vertical drift can reduce the impact caused by pre-sequence base addition and loss.

[0064] The method includes a core sequence set. The core sequence set is empty before clustering, and contains all original data sets after clustering. Therefore, the method constructs multiple tree structures in the front, middle and back of the core sequence. For each sequence in the core sequence set, the specified interval is added to the tree structure. The unclassified sequence is searched in the tree structure using the specified interval. If all nodes of the sequence are successfully searched, the sequence is matched to the core sequence. If the sequence cannot be searched, it is added to the core sequence set as a new sequence, and the specified interval is added to the tree structure. In addition, the node drift parameter is specified during searching, so that the same cluster sequence can be successfully searched for the error sequence.

[0065] Due to the characteristics of the tree structure, the time for searching the tree structure by the core sequence set is not affected. Obviously, the time complexity of the method is linear, that is, each sequence in the data set is executed once, and the processing time of the first sequence and the last sequence is theoretically equal. In addition, the algorithm can release the memory after reading the unclassified sequence, so the memory complexity is only related to the depth of the tree and the size of the original sequence, which greatly reduces the memory loss. In order to improve the practicability of the model, the unclassified sequence is allowed to be matched with the core sequence, and then a global comparison is performed to improve the accuracy of the core sequence. Then, the core sequence set after clustering can be directly outputted to simplify the reading of data. In addition, multiple processes are allowed to run to improve the execution speed of the algorithm. The specific method of multiple processes is that the sequences are divided into processes with the same initial base, and then an initial screening is performed.

[0066] As shown in Figure 1 , a flowchart of the method is given, including the following steps:

[0067] Step a: initialize parameters;

[0068] The method provided in the application allows customizing various model parameters, including but not limited to the number of tree structures, the length of the tree structure, the interval position selected by the tree structure, the longitudinal drift value, the threshold of the lateral drift value (i.e. if the lateral drift value is greater than the threshold, the search is abandoned), the number of processes (if greater than one, it is a multi-process mode, and the value can only be an exponential of 4, such as 4, 16, 64), the input file format (allowing input of labeled and unlabeled files, and allowing input of files in fasta, fastq, and txt formats), the output file format (allowing output of classification information of clusters, core sequence sets and other information, and if the input file is in a labeled format, it also allows output of accuracy, time consumption, coverage, etc. The core sequence set and the tree structure are constructed, and the constructed core sequence set and the tree structure will be empty sets, which will gradually increase only when the sequencing sequence enters the clustering operation. It should be noted that if it is a multi-process mode, the core sequence sets and the tree structures in different processes are not associated with each other. Only after all the sequences are clustered, the core sequence sets in different processes are combined.

[0069] Step c: determine whether the sequence to be tested is empty, if yes, go to step d; if no, go to step f;

[0070] Step d: output the cluster situation and the core sequence set;

[0071] Step e: end;

[0072] Step f: take out a sequence to be tested and search the tree structure;

[0073] Step g: determine whether a similar core sequence can be searched, if yes, go to step h; if no, go to step j;

[0074] Step h: classify the sequence to be tested and the similar core sequence into the same cluster;

[0075] Step i: optionally perform global alignment and correct core sequence errors;

[0076] Step j: add the sequence to be tested to the core sequence set and the tree structure, and go to step c.

[0077] It should be noted that for each sequence to be sequenced, it is searched against the tree structure, if a matching core sequence is successfully searched, it is divided into the cluster of the matched core sequence, if the global alignment function is turned on, the full sequence alignment is performed with the matched core sequence, and the positions with differences are marked, and then the core sequence is corrected. If the core sequence cannot be matched, it is added to the core sequence set as a new core sequence, and the specific interval is added to the tree structure. The specific algorithm details are as follows: first, compare the first interval with the first tree structure to see if they can be matched (i.e. the horizontal drift value generated by the matched index is less than the threshold value), if they can be matched, stop the subsequent tree structure search, if the global alignment function is turned on, at this time the sequence is globally aligned with the matched sequence (the algorithm and specific process of global alignment are not the focus of this patent, there are very mature global alignment algorithms at present, the program of this application allows to provide an interface, which can directly use the existing global alignment algorithm), after global alignment, the different base positions in the sequence will be marked, if a base position of a sequence in the core sequence set is frequently marked, it will be considered as an error base position and then corrected; if the first interval cannot be matched, then the end interval and the end tree structure are matched. If the end interval can be successfully matched, the optional global alignment and error correction function is performed as described above; if the end interval cannot be successfully matched, the middle interval is matched; when selecting the middle interval, the original specified interval is allowed to be shifted by the vertical drift value of base positions, and then multiple intervals are allowed to be selected, for example, the original middle interval is the 40th to the 60th base, if the vertical drift value is 2, the base interval of the final searched tree structure is [38, 58], [39, 59], [40, 60], [41, 61], [42, 62], after searching the series of base intervals of the tree structure, the interval with the smallest horizontal drift value is selected and matched with the matched sequence, if the horizontal drift value is less than the threshold value at this time, the matching is successful, and the optional global alignment and error correction function is performed; if it still cannot be matched, the sequencing sequence is regarded as a new cluster, which is added to the core sequence set, and its specified interval is added to the tree structure.

[0078] If it is a multi-process mode, after step a and before step c, step b is further included: the to-be-sequenced sequence is shunted. 1. If it is a multi-process mode, the data set is shunted to different processes. The specific shunting method is as follows: since the first interval error rate of the sequencing sequence is low, we will shunt according to the information of the first interval, for example, if the number of processes is 4, the first base of the sequencing sequence is shunted according to “C” “G” “T” “A”; if the number of processes is 16, the first two bases of the sequencing sequence are shunted. According to this shunting method, we can ensure that the sequences of the same cluster are divided into the same process, and ensure that a sequencing sequence cannot be divided into two different processes.

[0079] wherein, steps f-j specifically include:

[0080] Take out a sequence to be tested and search the tree structure to determine whether a similar core sequence can be searched. First, compare the first section interval with the first section tree structure to determine whether they can be matched. If they can be matched, stop the search of the subsequent tree structure, and classify the sequence to be tested and the similar core sequence as the same cluster. If the global alignment function is turned on, at this time, the sequence and the matched sequence are subjected to global alignment. After global alignment, the different base positions in the sequence will be marked. If a certain base position of a sequence in the core sequence set is frequently marked, it will be considered as an error base position and further corrected. If the first section interval cannot be matched, then the end section interval and the end tree structure are matched. If the end section interval can be successfully matched, stop the search of the subsequent tree structure, and classify the sequence to be tested and the similar core sequence as the same cluster. If the global alignment function is turned on, at this time, the sequence and the matched sequence are subjected to global alignment. After global alignment, the different base positions in the sequence will be marked. If a certain base position of a sequence in the core sequence set is frequently marked, it will be considered as an error base position and further corrected. If the end section also cannot be successfully matched, the middle section interval is matched. When selecting the middle section interval, the original specified interval is allowed to drift horizontally and vertically by a number of base positions, thereby allowing multiple intervals to be selected. After searching the base interval tree structure, the interval with the smallest horizontal drift value is selected and matched with the matched sequence. If the horizontal drift value at this time is less than the threshold value, the matching is successful. If the global alignment function is turned on, at this time, the sequence and the matched sequence are subjected to global alignment. After global alignment, the different base positions in the sequence will be marked. If a certain base position of a sequence in the core sequence set is frequently marked, it will be considered as an error base position and further corrected. If it still cannot be matched, the sequencing sequence is considered as a new cluster, which is added to the core sequence set, and its specified interval is added to the tree structure.

[0081] wherein, when all the sequences to be tested are completely clustered, a threshold value is allowed to be set. If a cluster contains fewer sequences, it is considered as a noise cluster, which is discarded from the core sequence set. In addition, if the multi-process mode is selected, the core sequence sets in different processes are merged, but different sequences in the sequence set are not merged. If the input file is a labeled data set, the time consumption, accuracy and correctness are calculated. Finally, the cluster results and the core sequence set are output.

[0082] In a preferred embodiment, the present application provides a nucleic acid sequence clustering device, comprising the following units:

[0083] A parameter initialization unit for initializing parameters;

[0084] A shunt unit for shunting the sequences to be tested;

[0085] a to-be-sequenced sequence judging unit, configured to judge whether the to-be-sequenced sequence is empty or not;

[0086] a result output unit, configured to output the cluster situation and the core sequence set;

[0087] a searching unit, configured to search a tree structure with a to-be-sequenced sequence;

[0088] a searching result judging unit, configured to judge whether a similar core sequence can be searched or not;

[0089] a judging result dividing unit, configured to divide the to-be-sequenced sequence and the similar core sequence into the same cluster, or add the to-be-sequenced sequence into the core sequence set and the tree structure;

[0090] a global alignment error correction unit, configured to perform global alignment and correct the core sequence error.

[0091] In a preferred embodiment, the present application provides a computer readable storage medium, comprising a program or instructions, when the program or instructions are executed on a computer, realizing any one of the nucleic acid sequence clustering methods as described above.

[0092] In a preferred embodiment, the present application provides a computer terminal, comprising a memory, and one or more processors in communication connection with the memory;

[0093] The memory stores instructions executable by the one or more processors, and the instructions are executed by the one or more processors to enable the one or more processors to realize the nucleic acid sequence clustering method as described above.

[0094] The application evaluates the performance of the method using real data and simulated data. Erlich and Zielinski et al. proposed a DNA encoding technology based on fountain code, which can recover several orders of magnitude more information than before. They synthesized 72000 DNA molecules with a length of 152. The DNA synthesis technology adopted the Twist company technology, and the sequencing adopted the MiseqV4 technology of the Illumina company. The application uses the real sequencing data set ERR181698 (a total of 14654644 sequences, belonging to 72000 original sequences) and ERR1817036 (a total of 34095791 sequences, belonging to 72000 original sequences) synthesized by them. In addition, starcode and DBSCAN are selected as benchmark algorithms, where starcode is the fastest DNA clustering algorithm in the currently disclosed algorithm, which mainly calculates the Levenshtein distance by editing the matrix. DBSCAN is a traditional clustering algorithm with low complexity, which is a density-based clustering algorithm with less than quadratic computational complexity. The test environment includes a home computer, a cloud server, and a supercomputer. In order to obtain a high-credibility labeled data set, we use biological software such as pear and bowtie to compare the original data set and the original set.

[0095]

[0096] Table 1: Time consumption comparison of the method (unit: seconds)

[0097] The ERR181698 data set is extracted as a benchmark data set of different orders of magnitude, to more quantitatively show the clustering effect. Table 1 shows the time consumption comparison between different algorithms. From the table, it can be seen that the clustering speed of the method of the application is higher than that of starcode and DBSCAN, and multi-core running will significantly improve the clustering speed of the algorithm. In addition, Table 2 shows the accuracy of the method on the real data set, which can be seen that the method has very high accuracy on the real data set.

[0098]

[0099] Table 2: Accuracy of the method

[0100] In order to show the clustering effect of the method on a large amount of data set, a set of 10 billion DNA sequencing sequences (sequencing depth is 1000, error rate is 0.04%) is simulated, which is the largest simulated data set in the field of DNA storage known at present. The experiment is carried out in a supercomputer server in multiple threads, and the experimental results show that the method can cluster 1 billion data in about 4 hours under multi-process, and the clustering accuracy is 99.99%. It is proved that the method can still cluster in a short time when facing large-scale data set.

[0101] It should be noted that the technical solutions not described in detail in the present application adopt known technologies.

[0102] The above only describes the preferred embodiments of the present application, and it should be pointed out that for ordinary skilled in the art, without departing from the principles of the present application, a number of improvements and refinements can be made, and these improvements and refinements should be considered as the protection scope of the present application.

Claims

1. A nucleic acid sequence clustering method, characterized in that, Includes the following steps: Step a: Initialize parameters; Step c: Determine if the sequence to be tested is empty. If yes, proceed to step d; otherwise, proceed to step f. Step d: Output cluster information and core sequence set; Step e: End; Step f: Extract a sequence to be tested and perform a search on the tree structure; Step g: Determine if a similar core sequence can be retrieved. If yes, proceed to step h; otherwise, proceed to step j. Step h: Group the sequence to be tested with similar core sequences into the same cluster; Step i: Perform global alignment and correct core sequence errors; Step j: Add the sequence to be tested to the core sequence set and the tree structure, then skip to step c; Step fj specifically includes: Step 1: Take out a sequence to be tested and search it against the tree structure to determine whether similar core sequences can be found; Step 1.a: First, compare whether the first segment interval can match the first segment tree structure. If they can match, stop the search of subsequent tree structures and classify the test sequence and similar core sequences into the same cluster. Step 1.b: If the first interval cannot be matched, then match the last interval with the last tree structure. If the last interval is successfully matched, stop the search of the subsequent tree structure and classify the test sequence and similar core sequences into the same cluster. Step 1.c: If the end also fails to match, then the middle interval is matched. When selecting the middle interval, the horizontal and vertical drift of several base positions will be allowed on the basis of the original specified interval, thus allowing the selection of multiple intervals. After the base interval search tree structure is completed, the interval with the smallest horizontal drift value and the matching sequence will be selected. If the horizontal drift value is less than the threshold at this time, the match is successful, and the sequence to be tested and similar core sequences are classified into the same cluster. Step 2: If the global alignment function is enabled, the sequence will be compared with the matched sequence. After the alignment, different base positions in the sequence will be marked. If a certain base position of a sequence in the core sequence set is frequently marked, it will be regarded as an erroneous base position and then corrected. Step 3: If a match still cannot be found, treat the sequencing sequence as a new cluster, add it to the core sequence set, and add its specified interval to the tree structure.

2. The nucleic acid sequence clustering method according to claim 1, characterized in that, In step a, the initialization parameters include the number of initial tree structures, the length of the tree structure, the selected interval position of the tree structure, the vertical drift value, the horizontal drift value threshold, the number of processes, the input file format, the output file format, the core sequence set, and the tree structure.

3. The nucleic acid sequence clustering method according to claim 1, characterized in that, If it is a multi-process mode, after step a and before step c, it also includes step b: splitting the test sequence.

4. The nucleic acid sequence clustering method according to claim 1, characterized in that, Between steps c and d, there is also step k, which sets a threshold. If a cluster contains few sequences, it is considered a noisy cluster and is discarded from the core sequence set. In addition, if the multi-process mode is selected, the core sequence sets of different processes are merged, but different sequences in the sequence set are not merged. If the input file is a labeled dataset, the time consumption, accuracy, and correctness are calculated. Finally, the cluster results and the core sequence set are output.

5. A nucleic acid sequence clustering device, characterized in that, Includes the following units: The parameter initialization unit is used to initialize parameters. The splitting unit is used to split the sequence under test; The test sequence judgment unit is used to determine whether the test sequence is empty; The result output unit is used to output cluster information and the core sequence set; The retrieval unit is used to retrieve a sequence to be tested and perform a search on the tree structure. The retrieval result judgment unit is used to determine whether similar core sequences can be retrieved; The judgment result is divided into units, which are used to classify the test sequence and similar core sequences into the same cluster, or to add the test sequence to the core sequence set and tree structure; The global alignment and error correction unit is used to perform global alignment and correct core sequence errors; The retrieval result judgment unit, the judgment result partitioning unit, and the global comparison and error correction unit are specifically used to perform the following steps: Step 1: Take out a sequence to be tested and search it against the tree structure to determine whether similar core sequences can be found; Step 1.a: First, compare whether the first segment interval can match the first segment tree structure. If they can match, stop the search of subsequent tree structures and classify the test sequence and similar core sequences into the same cluster. Step 1.b: If the first interval cannot be matched, then match the last interval with the last tree structure. If the last interval is successfully matched, stop the search of the subsequent tree structure and classify the test sequence and similar core sequences into the same cluster. Step 1.c: If the end also fails to match, then the middle interval is matched. When selecting the middle interval, the horizontal and vertical drift of several base positions will be allowed on the basis of the original specified interval, thus allowing the selection of multiple intervals. After the base interval search tree structure is completed, the interval with the smallest horizontal drift value and the matching sequence will be selected. If the horizontal drift value is less than the threshold at this time, the match is successful, and the sequence to be tested and similar core sequences are classified into the same cluster. Step 2: If the global alignment function is enabled, the sequence will be compared with the matched sequence. After the alignment, different base positions in the sequence will be marked. If a certain base position of a sequence in the core sequence set is frequently marked, it will be regarded as an erroneous base position and then corrected. Step 3: If a match still cannot be found, treat the sequencing sequence as a new cluster, add it to the core sequence set, and add its specified interval to the tree structure.

6. A computer-readable storage medium, characterized in that, Includes a program or instructions that, when run on a computer, implement the nucleic acid sequence clustering method as described in any one of claims 1-4.

7. A computer terminal, characterized in that, Includes a memory, and one or more processors communicatively connected to the memory; The memory stores instructions that can be executed by the one or more processors to enable the one or more processors to implement the nucleic acid sequence clustering method as described in any one of claims 1-4.

Images